Climate pathways behind phytoplankton-induced atmospheric warming

. We investigate the ways in which marine biologically mediated heating increases the surface atmospheric temperature. While the effects of phytoplankton light absorption on the ocean have gained attention over the past years, the impact of this biogeophysical mechanism on the atmosphere is still unclear. Phytoplankton light absorption warms the surface of the ocean, which in turn affects the air–sea heat and CO 2 exchanges. However, the contribution of air– sea heat versus CO 2 ﬂuxes in the phytoplankton-induced atmospheric warming has not been yet determined. Different so-called climate pathways are involved. We distinguish heat exchange, CO 2 exchange, dissolved CO 2 , solubility of CO 2 and sea-ice-covered area. To shed more light on this subject, we employ the EcoGEnIE Earth system model that includes a new light penetration scheme and isolate the effects of individual ﬂuxes. Our results indicate that phytoplankton-induced changes in air–sea CO 2 exchange warm the atmosphere by 0.71 ◦ C due to higher greenhouse gas concentrations. The phytoplankton-induced changes in air–sea heat exchange cool the atmosphere by 0.02 ◦ C due to a larger amount of outgoing longwave radiation. Overall, the enhanced air–sea CO 2 exchange due to phytoplankton light absorption is the main driver in the biologically induced atmospheric heating.


Introduction
Previous studies have shown that marine biota can modify the light penetration in the ocean with consequences on the atmospheric temperature and on the climate system (Shell et al., 2003;Wetzel et al., 2006;Gnanadesikan and Anderson, 2009).Using an Earth system model (ESM) of intermediate complexity, we identify and compare the climate pathways behind the changes in atmospheric temperature due to phytoplankton light absorption.
Marine biota and phytoplankton play a major role in the absorption of light and therefore in the vertical distribution of heat in the upper layers of the ocean (Kowalczuk et al., 2019).Indeed, observational evidence supports the hypothesis that chlorophyll increases the upper ocean heat uptake.For instance, satellite observations show that phytoplankton blooms can cause an increase in sea surface temperature (SST) of 1.5 • C (Kahru et al., 1993).Furthermore, previous remote sensing data indicate an increase in local SST of 4.5 • C on a 4 d timescale due to the presence of phytoplankton blooms (Capone et al., 1998).Recent high-resolution in situ observations in the Indo-West Pacific Ocean highlight large anomalies of temperature of 0.95 • C in the uppermost skin layer of the ocean when large phytoplankton blooms appear (Wurl et al., 2018).However, all these observations are either on a short timescale or in a geographically limited area.
To study the larger-scale impact of phytoplankton light ab-sorption and its relative magnitude, Earth system models are employed.
Models of differing complexity are used to study the effect of phytoplankton light absorption.For instance, using oceanonly (Anderson et al., 2007) or general circulation models, several studies focusing on the tropical Pacific Ocean (Murtugudde et al., 2002;Lengaigne et al., 2007;Löptien et al., 2009) or on the Arctic Ocean (Lengaigne et al., 2009) report an increase in SST between 0.5-2 • C due to phytoplankton light absorption.A warming of the ocean surface induced by marine biota has consequences on the overall climate system.For instance, Patara et al. (2012) find that an increase in SST due to phytoplankton light absorption increases the atmospheric humidity content thereby increasing the greenhouse effect and the atmospheric temperature locally by up to 0.5 • C. Furthermore, phytoplankton can amplify locally the seasonal temperature of the lowest atmospheric layer by 1 • C, changing the Walker and Hadley circulation (Shell et al., 2003).
It is therefore known that phytoplankton light absorption plays a non-negligible role in the atmospheric temperature, but which climate pathway is the most important behind this warming is still unclear.Phytoplankton light absorption affects the surface atmospheric temperature via two climate pathways.First, various modeling studies suggest that biologically induced surface water heating can increase the airsea heat exchange (Capone et al., 1998;Oschlies, 2004;Wetzel et al., 2006) with consequences on the formation of tropical storms and monsoons in the Arabian Sea (Sathyendranath et al., 1991).Second, the solubility of gases and thus also the air-sea CO 2 exchange is affected by phytoplankton light absorption.For instance, Manizza et al. (2008) study the impact of this biogeophysical mechanism on the air-sea flux of CO 2 and find that phytoplankton light absorption has a small outgassing effect on a global scale with high regional fluctuations.
However, none of these studies have analyzed, disentangled and compared the changes in both air-sea heat and CO 2 exchange due to phytoplankton light absorption.To shed light on the biologically induced atmospheric warming, we use a recent Earth system model called EcoGEnIE (Ward et al., 2018).In an earlier study, we implemented phytoplankton light absorption in this model (Asselot et al., 2021).We consider two different biologically induced changes: changes in air-sea heat and changes in air-sea CO 2 exchange (Fig. 1).The air-sea CO 2 exchange can be influenced by the dissolved oceanic CO 2 in three different ways: (1) through changes in the biogeochemical pumps as a result of phytoplankton light absorption -for instance, Manizza et al. (2008) have shown that changes in oceanic circulation due to phytoplankton light absorption enhance the vertical supply of nutrients, increasing the relative abundance of calcifiers, and as a consequence, the primary production and the export production of organic matter increase; (2) through a decrease in CO 2 solubility due to higher SST, increasing the atmospheric CO 2 concentrations and the greenhouse gas effect; and (3) through a decrease in sea-ice formation, because sea ice acts as an ocean cap that blocks gas exchanges.To achieve the disentangling of the specific climate pathways, we turn them on and off by prescribing values in our ESM in order to isolate their impact on the climate system.

Model description
Our motivation is to study the interactions between the marine ecosystem, the biogeochemistry, the biogeophysics and the climate system.These interactions are computationally expensive in high-resolution models; therefore we used an Earth system model of intermediate complexity (Claussen et al., 2002).The Earth system model employed is the carbon-centric Grid-Enabled Integrated Earth system model (cGEnIE) (Lenton et al., 2007) composed of several modules describing the dynamics of the Earth system (Fig. 2).This model has been previously calibrated and compared to observations (Edwards and Marsh, 2005;Lenton et al., 2006;Ridgwell et al., 2007;Marsh et al., 2011).Moreover, this model is widely used to study past climate and changes in the carbon cycle over geological times (Ödalen et al., 2018;Greene et al., 2019;Adloff et al., 2020), past mass extinctions (Alvarez et al., 2019) and biogeochemistry processes (Meyer et al., 2016).Additionally, cGEnIE has been employed to assess the sensitivity of atmospheric CO 2 to biogeochemical pumps, ocean circulation and climate feedbacks in the Southern Ocean (Cameron et al., 2005).A new ecosystem component (ECOGEM) is associated with cGEnIE to form the recent EcoGEnIE model (Ward et al., 2018).EcoGEnIE is used to determine the link between the marine plankton ecosystem and various past climate scenarios (Wilson et al., 2018) with a focus on phosphorus inventory (Reinhard et al., 2020).For our study, the model combines different components including ocean hydrodynamics, atmosphere, sea ice, ocean biogeochemistry and marine ecosystem.We do not consider a terrestrial component meaning that the land surface is essentially passive.We use the same configuration as described in detail by Asselot et al. (2021), and the following description only refers to our specific model setup.

The physical components
The physics of the model contains a frictional-geostrophic ocean circulation (GOLDSTEIN), coupled to a 2D energymoisture balance model of the atmosphere (EMBM) and a thermodynamic sea-ice model (GOLDSTEINSEAICE) (Edwards and Marsh, 2005;Marsh et al., 2011).Heat and moisture are exchanged between the three components and act as a coupling strategy.The oceanic component calculates the horizontal and vertical redistribution of heat, salinity and biogeochemical elements via advection, convection and mixing.The ocean module is configured on a 36 × 36 horizontal grid.The horizontal grid is 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 grid has been used for previous biogeochemical simulations (Cameron et al., 2005;Colbourn, 2011).We consider 32 vertical oceanic layers increasing logarithmically from 29.38 m for the surface layer to 456.56 m for the deepest layer.This vertical resolution is already used to study the relative importance of biogeophysical and biogeochemical mechanisms on the climate system (Asselot et al., 2021).The atmospheric component is based closely on the UVic Earth system model (Weaver et al., 2001).The prognostic variables are atmo-spheric temperature and specific humidity.Precipitation removes instantaneously all moisture corresponding to an excess above a relative humidity threshold.The wind stress is prescribed and identical between all simulations; the temperature cannot affect the wind stress.The sea-ice component solves the equation for part of the ocean covered by sea ice.The prognostic variables are ice thickness and ice areal fraction.The transport of sea ice includes sources and sinks of these variables.The growth or decay of sea ice depends on the net heat flux into the ice.The dynamics in this module consist of advection by currents and diffusion.Sea ice does not limit the penetration of photosynthetically available radiation in the ocean.The biogeochemical module (BIOGEM) represents the transformation and spatial redistribution of biogeochemical tracers (Ridgwell et al., 2007).The state variables are inorganic and organic matter.The biological uptake is represented by an implicit biological community: nutrients are directly converted into organic matter via an uptake rate.The biological uptake is limited by light, temperature and nutrient availability.Organic matter is partitioned into dissolved and particulate phases (DOM and POM).For this study, BIO-GEM does not consider a temperature dependency on iron solubility and iron bioavailability.Our model setup includes iron (Fe) and phosphate (PO 4 ) as limiting nutrients.Similar to Asselot et al. (2021), we do not consider nitrate here.Moreover, the surface production is redistributed in the water column as a depth-dependent flux.To achieve this, the surface export is divided between refractory organic matter remineralized close to the seafloor and labile organic matter mostly remineralized in the upper water column (Ridgwell et al., 2007).Furthermore, because we do not consider a sediment component, all organic matter reaching the sea floor is instantaneously remineralized.Calcium carbonate (CaCO 3 ) is represented in the model, and its dissolution below the surface follows the same profile as the remineralization of POM.Additionally, CaCO 3 production is parameterized as a fixed ratio to POC production and so is not independent in our model setup.Recent studies have implemented and calibrated a temperature-dependent remineralization in the model (Crichton et al., 2021;Armstrong McKay et al., 2021), but this parameterization is not included in our model setup.Furthermore, BIOGEM calculates the air-sea CO 2 and O 2 exchange.The value of atmospheric CO 2 predicted by BIO-GEM is used as input for the radiative scheme of the atmospheric component, thus providing climate feedback.

Ecosystem component
The marine ecosystem component (ECOGEM) represents the marine plankton community and associated interactions within the ecosystem (Ward et al., 2018).The biological uptake in ECOGEM replaces the BIOGEM uptake calculation and is limited by light, temperature and nutrient availability.Plankton biomass and organic matter are subject to processes such as resource competition and grazing before being passed to DOM and POM.Several ecophysiological parameters are size-dependent such as maximum nutrient uptake rate, cell carbon quotas, grazing, and partitioning between DOM and POM.Additionally, the nutrient uptake, photosynthesis and predation are temperature-dependent.The ecosystem is divided into different plankton functional types (PFTs) with specific traits.Each PFT is sub-divided into size classes with specific size-dependent traits.Here, we consider only two PFTs: phytoplankton and zooplankton (Appendix A1).
Phytoplankton is characterized by nutrient uptake and pho-tosynthesis, whereas zooplankton is characterized by predation traits.Zooplankton grazing depends on the concentration of prey biomass, with predominant grazing on prey that are 10 times smaller than themselves.Each population is associated with biomass state variables for carbon, phosphate, iron and chlorophyll.The production of dead organic matter is a function of mortality and messy feeding, with partitioning between non-sinking dissolved and sinking particulate organic matter.Finally, plankton mortality is reduced at very low biomass such that plankton cannot become extinct.

Light absorption in the ocean
In the previous model version (Ward et al., 2018), light was only absorbed by phytoplankton.In the model version of Asselot et al. ( 2021), a new light scheme is implemented where the absorbed light by phytoplankton is converted into heat and is able to affect the oceanic temperature.Furthermore, light absorption takes place throughout the water column and is not restricted to the first oceanic layer anymore.The light absorption scheme is a coupling between Eqs. ( 1) and (2).For simplicity, in our model configuration, the incoming shortwave radiation does not vary seasonally.We look at long-term changes in the climate system; therefore the absence of a seasonal cycle does not affect the overall qualitative and main findings.The presence of organic and inorganic particles as well as dissolved molecules restrains the light penetration in the ocean (Ward et al., 2018).The vertical light attenuation scheme is given by Eq. (1): where I (z) is the irradiation of the full solar spectrum at depth z, I 0 is the irradiation at the surface of the ocean, k w is light absorption by clear water and inorganic particles (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 taken from Ward et al. (2018).The parameter I 0 is negative in the model, because it is a downward flux from the sun to the surface of the ocean.We allow primary production and light to penetrate until the sixth layer of the model (221.84m deep), which is the lower limit of the euphotic zone (Tett, 1990).In our model setup, maximum absorption occurs in the upper oceanic layer and minimum absorption occurs in the sixth layer.
Phytoplankton changes the optical properties of the ocean (Sonntag and Hense, 2011) through phytoplankton light absorption.Therefore it can cause a radiative heating and change the oceanic temperature.We implemented phytoplankton light absorption into the model following Hense (2007) and Patara et al. (2012).The scheme is given by Biogeosciences , 19, 223-239, 2022 https://doi.org/10.5194/bg-19-223-2022 Eq. ( 2): (2) ∂T /∂t denotes the temperature changes, c p is the specific heat capacity of water, ρ is the ocean density and I is the solar radiation incident at depth z.Part of the light absorbed is used by phytoplankton for photosynthesis and part leads to heating of the water.

Air-sea heat exchange
We detail here the total heat flux from the ocean and sea ice going into the atmosphere.The vertically integrated atmospheric heat equation (Eq. 3) is given by Weaver et al. (2001) and Marsh et al. (2011): Q ta corresponds to the total heat flux into the atmosphere, Q SW is the net shortwave radiation corresponding to the solar irradiance received from the sun and reflected by the planet's albedo, C A is a heat absorption coefficient (0.3 over the ocean; Marsh et al., 2011), Q LH is the latent heat flux corresponding to phase change of a thermodynamic system, Q SH is the sensible heat flux corresponding to temperature change of a thermodynamic system, Q LW is the net (upward minus downward) re-emitted longwave radiation corresponding to infrared energy coming from the planet and Q PLW is the outgoing planetary longwave radiation.The atmosphere loses heat through net longwave radiation, dominated by the outgoing longwave radiation; thus the total longwave heat flux (Q LW + Q PLW ) is negative in the model.Furthermore, evaporative cooling of the ocean leads to a latent heat release in the atmosphere upon condensation and precipitation.Evaporated water vapor may be transported away from an oceanic source, to condense and precipitate elsewhere.

Air-sea CO 2 exchange
The atmospheric temperature depends on the atmospheric CO 2 concentration which is affected by the transfer of CO 2 between the ocean and the atmosphere.The flux of CO 2 across the atmosphere-ocean interface (Eq. 4) is given by Ridgwell et al. (2007): (4) F CO 2 is the air-sea CO 2 flux, k corresponds to the gas transfer velocity, ρ is the ocean density, C w is the concentration of dissolved gas in the surface ocean, α is the solubility coefficient calculated from Wanninkhof (1992) and depends on the sea surface temperature and salinity, C a is the concentration of gas in the atmosphere, and A is the fraction of the ocean covered by sea ice.
Phytoplankton light absorption warms the surface of the ocean and thus reduces CO 2 solubility and sea-ice fraction.The flux of CO 2 is therefore affected via the parameters C w , α and A. To study the flux precisely we either prescribe these parameters in the air-sea CO 2 exchange calculation or let them evolve freely (see below).

Model setup and simulations
During this study, we are primarily interested in the relative differences between our selected simulations.We focus on the relative impact of phytoplankton light absorption on different climate pathways rather than on the changes in the climate state.We try to simulate realistic mean climate systems but the absolute values of the climate quantities are less relevant due to the limitations of such a model of intermediate complexity.
For a realistic nutrient distribution in the ocean, we performed a BIOGEM spin-up for 10 000 years.During the spin-up the atmospheric CO 2 concentration is fixed to 278 ppm.The simulations restart for 1000 years after the spin-up with ECOGEM, meaning that all simulations consider marine biota.Due to the single-layer atmospheric component, the non-seasonality and the non-representation of the land dynamics, running the model for 1000 years is sufficient to achieve a steady state.The results represent the annual mean of the last year of the simulations, when the model is in a steady state.The present-day continental configuration, model setup, grid resolution and ecosystem community are identical as in Asselot et al. (2021).Our model setup has a reasonable match to observations, and further details can be found in Asselot et al. (2021).For simplification, only one phytoplankton and one zooplankton species are included in the model setup (Appendix A1 and B1).Repeating our main simulations with multiple size classes leads to small differences compared to the simulations with one size class (Appendix B1).Unlike our previous study, prescribing SST is necessary for some simulations here.Due to technical issues, we cannot prescribe the seasonal cycle of SST but only the annually averaged SST.As a consequence, the seasonality is removed from our model setup.The absence of the seasonal cycle is not an issue for this study because we look at the importance of each climate pathway rather than focusing on the quantitative changes of the climate system.
The carbon cycle is closed in our simulations, meaning that there is no input of carbon through volcanic or anthropogenic activities.Only the size of the carbon reservoirs can vary.If not stated otherwise, the concentration of atmospheric CO 2 evolves freely in the simulations.All simulations are forced with the same constant flux of dissolved iron into the ocean surface (Mahowald et al., 2006).To ensure that the model is suitable for our study, we conducted two sensitivity analyses.First, we analyzed the sensitivity of the climate variables by conducting two simulations with differhttps://doi.org/10.5194/bg-19-223-2022 Biogeosciences, 19, 223-239, 2022 ent atmospheric CO 2 concentrations (Appendix C1).Second, we ensure that the heat and CO 2 interaction is negligible (Appendix D1).To study the effect of phytoplankton light absorption on the atmospheric temperature we perform seven different simulations, all including the ECOGEM component (Fig. 3): -Bio is the reference run and is the only simulation that does not include phytoplankton light absorption (k Chl = 0 in Eq. 1).In this simulation, all the climate pathways evolve freely.
-BioLA is the same as the reference run but phytoplankton light absorption is implemented.In this simulation, all the climate pathways evolve freely.
-HEAT is the same as the second one except that we prescribe the atmospheric CO 2 concentration only for the atmospheric temperature calculation.For a comparison with the reference run, the prescribed atmospheric CO 2 concentration from Bio is used (169 ppm).The effect of CO 2 on atmospheric temperature is fixed but the airsea heat fluxes evolve freely.This simulation analyses the effect of phytoplankton-induced changes of air-sea heat fluxes on the atmospheric temperature.
-CARB is the simulation with an uncoupled oceanatmosphere setup.The atmospheric component is forced with the heat fluxes from Bio and the atmospheric CO 2 concentration is prescribed with the value of BioLA.This simulation determines the effect of phytoplankton-induced changes of atmospheric CO 2 concentration on the atmospheric temperature.Please note that CARB is well suited for studying the atmospheric properties but not to examine ocean dynamics.
-HCorg is the simulation where we only allow the biogeochemical pumps (soft-tissue pump and carbonate pump) to affect the dissolved CO 2 .The solubility of CO 2 (α in Eq. 4) and sea-ice extent (A in Eq. 4) parameters are prescribed using the respective values from Bio.
The CO 2 solubility is fixed by prescribing the SST only for this calculation.In HCorg air-sea heat exchange and the biogeochemical pumps parameter (C w in Eq. 4) evolve freely.
-HCorgSI is the simulation where the biogeochemical pumps and sea-ice extent affect dissolved CO 2 .The CO 2 solubility (α in Eq. 4) is prescribed using the value of Bio.In HCorgSI the air-sea heat exchange, the biogeochemical pumps (C w in Eq. 4) and sea-ice extent (A in Eq. 4) parameters evolve freely.
-HCorgSol is the simulation where the biogeochemical pumps and the solubility pump affect dissolved oceanic CO 2 .The sea-ice extent parameter (A in Eq. 4) is prescribed using the value of Bio.In HCorgSol the air-sea heat exchange, the biogeochemical pumps (C w in Eq. 4) and the CO 2 solubility (α in Eq. 4) parameters evolve freely.
4 Global response of the climate system In this section we present the results of the simulations on a global scale.We do not consider local patterns because we removed any seasonal cycle in our model setup.Moreover, the horizontal grid resolution is low and marine biota cannot move between grid cells; thus even if seasonality was included, key regional patterns will not be resolved.First, we focus on the chlorophyll biomass and sea surface temperature, because phytoplankton light absorption has a direct effect on these climate variables (Oschlies, 2004;Lengaigne et al., 2007;Paulsen et al., 2018).Second, these changes in oceanic properties affect the carbon cycle (Manizza et al., 2008;Asselot et al., 2021); therefore we study the changes in atmospheric CO 2 concentration.Third, phytoplankton light absorption alters the atmospheric properties (Patara et al., 2012); thus we analyze the changes in radiative heat fluxes, humidity and evaporation.Finally, the response of the surface atmospheric temperature is analyzed.

Chlorophyll biomass and sea surface temperature
Our results indicate differences of SST and chlorophyll biomass, depending on the climate pathways included in our model setup (Table 1).Due to the uncoupled oceanatmosphere setup in CARB, ocean dynamics are not presented in this section.The reference run Bio has the lowest chlorophyll biomass and a low SST while the simulation BioLA has the highest chlorophyll biomass and SST.The increase in chlorophyll biomass is due to two different mechanisms: first, phytoplankton light absorption leads to a higher surface production, enhancing the remineralization at the surface of the ocean as shown by Asselot et al. (2021).Second, phytoplankton light absorption enhances the upward vertical velocity in the upwelling regions.As a result of these two mechanisms, the surface nutrient concentrations increase, explaining the higher chlorophyll biomass.The increase in surface chlorophyll biomass due to phytoplankton light absorption between BioLA and Bio is 0.012 mgChl m −3 , in line with a previous estimate of 0.014 mgChl m −3 (Asselot et al., 2021).The higher chlorophyll biomass is, however, limited by the increase in zooplankton biomass applying a top-down control.Via the effect of phytoplankton light absorption, a higher surface chlorophyll biomass leads to an increase in SST.The global difference of SST between BioLA and Bio is 0.08 • C.This value is lower than previous model estimates that show a global SST increase of 0.45-1 • C due to phytoplankton light absorption (Murtugudde et al., 2002;Löptien et al., 2009;Asselot et al., 2021).This underestimation of the biologically induced SST increase is due to the non-  seasonal radiative forcing of the model, decreasing the global heat budget (Appendix E1).
In HEAT, the chlorophyll biomass is higher while the SST is lower compared to Bio (Table 1).This is rather counterintuitive and is due to changes in oceanic circulation between these two simulations.For instance, the maximum Atlantic meridional overturning circulation (AMOC) is 8.6 Sv in HEAT while it is 7.6 Sv in Bio.In HEAT, the SST is lower and the sea-ice cover is slightly higher (Appendix F1) compared to Bio, leading to more deep water formation in polar latitudes.As a result, the AMOC is enhanced in HEAT.The stronger oceanic circulation in HEAT leads to an enhanced nutrient redistribution, thus increasing the surface nutrient concentrations.For instance, the surface PO 4 concentration is ∼ 0.21 µmol kg −1 in HEAT while it is ∼ 0.19 µmol kg −1 in Bio.The higher surface PO 4 concentration in HEAT explains the higher chlorophyll biomass in this simulation compared to Bio.The enhanced oceanic circulation in HEAT compared to Bio also leads to a stronger redistribution of heat along the water column, explaining the surface cooling and the warming of the deep ocean in HEAT.Our results indicate that the bottom water temperature in HEAT is 3.57 • C while it is 3.09 • C in Bio.
The surface chlorophyll biomass in the simulations HCorg, HCorgSI and HCorgSol are higher than the surface chlorophyll biomass in Bio due to the higher surface PO 4 concentrations.In Bio the surface PO 4 concen-tration is 0.18 µmol kg −1 while in the simulations HCorg, HCorgSI and HCorgSol the surface PO 4 concentrations are > 0.21 µmol kg −1 .The higher surface PO 4 concentrations are due to enhanced remineralization at the ocean surface and enhanced upward vertical velocities in the upwelling regions.Due to the effect of phytoplankton light absorption, the higher surface chlorophyll biomasses in HCorg, HCorgSI and HCorgSol lead to higher SSTs compared to Bio.Only the sea-ice extent differs between the simulations HCorg and HCorgSI, but their chlorophyll biomass and SSTs are identical.This result evidences a lack of sea-ice influence on these climate variables and thus on the heat fluxes.In addition, the chlorophyll biomass and SST are higher in HCorg compared to HCorgSol, indicating that the solubility factor has a negative effect on these climate variables.Between these two simulations, the only difference is the CO 2 -solubility factor that can evolve freely in HCorgSol.In the simulation HCorg, the SST for the calculation of the CO 2 solubility is prescribed using the value of Bio, which is the lowest value.Considering the physical and chemical properties of the ocean, a low SST increases the solubility of CO 2 (Wanninkhof, 1992).Therefore, the CO 2 solubility is reduced in HCorgSol compared to HCorg, due to the higher SST in HCorgSol.For instance, our results indicate that on a global scale, the surface oceanic CO 2 concentration is 27.200 µmol kg −1 in HCorgSol while it is 27.213 µmol kg −1 in HCorg.Via the nutrient ratios, these changes in carbon cycle between the https://doi.org/10.

Atmospheric properties
The oceanic properties differ between the simulations; thus we expect differences in the atmospheric properties.We compare the atmospheric CO 2 concentration, the heat fluxes, the evaporation, the specific humidity and finally the surface atmospheric temperature between the simulations.

Atmospheric CO 2 concentration
The atmospheric CO 2 concentrations for the simulations is low compared to the pre-industrial level (Fig. 4).This is due to our new model setup that allows primary production until the sixth oceanic layer, meaning that more carbon is stored in the deep ocean, reducing the atmospheric CO 2 concentration (see Asselot et al., 2021).In all the simulations considering phytoplankton light absorption, the atmospheric CO 2 concentration is higher than in the reference run.In a previous study, we evidence that the higher atmospheric CO 2 concentration is mainly due to a decrease in CO 2 solubility via the higher SST while the enhanced remineralization of organic matter and the dissolution of CaCO 3 affect the atmospheric CO 2 concentration slightly (Asselot et al., 2021).The atmospheric CO 2 concentration is the lowest in Bio while it is the highest in BioLA, with a difference of 9 ppm.This value is lower than a previous estimate that indicates an increase in atmospheric CO 2 concentration of 18 ppm (Asselot et al., 2021).This lower estimate is due to the non-seasonal cycle forcing, neglecting the seasonal variations of air-sea CO 2 exchanges.
In HEAT, the atmospheric CO 2 concentration is prescribed only for the atmospheric temperature calculation.Therefore the atmospheric CO 2 concentration can vary due to changes in dissolved oceanic CO 2 , sea-ice extent and CO 2 solubility, affecting the other climate variables.The atmospheric CO 2 concentration in HEAT is slightly higher than in Bio.This is due to the larger chlorophyll biomass in HEAT than in Bio (Table 1), indicating a higher production and thus more remineralization in the ocean.During the remineralization process, CO 2 is produced; thus the higher remineralization in HEAT increases the dissolved CO 2 concentration.On a global scale, our results indicate that the surface dissolved oceanic CO 2 is about 6.354 mol kg −1 in HEAT, while it is 6.302 mol kg −1 in Bio.The larger dissolved oceanic CO 2 concentration in HEAT increases the air-sea CO 2 flux and in turn the atmospheric CO 2 concentration (see Eq. 4).
The atmospheric CO 2 concentration in CARB is similar to the one in BioLA, because we prescribed the value against the one in BioLA.
The simulations HCorg, HCorgSI and HCorgSol have a higher atmospheric CO 2 concentration than in Bio.This is again not surprising because these simulations consider phytoplankton light absorption, which increases the atmospheric CO 2 concentration as shown by Asselot et al. (2021).The atmospheric CO 2 concentration between HCorg and HCorgSI is similar even if their sea-ice extent and sea-ice thickness differs (Appendix F1).The changes in sea ice do not have an effect on the atmospheric CO 2 concentration.The slightly higher atmospheric CO 2 concentration in HCorgSol compared to HCorg is due to changes in CO 2 solubility: as described above, the CO 2 solubility is lower in HCorgSol compared to HCorg.As a consequence, the air-sea CO 2 flux is higher in HCorgSol compared to HCorg, leading to a slightly higher atmospheric CO 2 concentration in HCorgSol.

Heat fluxes
The air-sea heat flux is divided into the net shortwave radiation, the net re-emitted longwave radiation, the sensible heat flux and the latent heat flux (Fig. 5).The air-sea heat fluxes represent the total heat fluxes from the ocean and sea ice, going into the atmosphere.The simulations HCorg and HCorgSI have exactly the same heat fluxes.The only difference between these two simulations is the different seaice extent for the calculation of the air-sea CO 2 flux.This change in air-sea CO 2 flux does not alter the air-sea heat flux, explaining the identical radiative heat fluxes between HCorg and HCorgSI.Furthermore, the simulations BioLA and HCorgSol have the same heat fluxes, and the only difference is also the sea-ice extent.As detailed previously, the changes in sea-ice extent do not affect the heat fluxes, explaining the identical radiative fluxes between BioLA and HCorgSol.Finally, the heat fluxes between CARB and Bio are identical because we prescribed the heat fluxes in CARB with the values of Bio.
The net shortwave heat flux is divided into two fluxes: the incoming shortwave radiation from the sun entering the  atmosphere and the outgoing reflected shortwave radiation leaving the atmosphere.Figure 5a shows that the net shortwave heat flux is identical for all the simulations and is positive.The positive values indicate that net shortwave heat flux is dominated by the flux entering the system, the incoming radiation.The incoming shortwave radiation from the sun is always identical between simulations.Therefore identical net shortwave heat flux implies that the outgoing reflected shortwave radiation is also the same between simulations due to the treatment of shortwave radiation in the model (Weaver et al., 2001).The net longwave heat flux is negative for all simulations, indicating that this flux is dominated by the upward longwave radiation leaving the atmosphere (Fig. 5b).A more negative value of net longwave heat flux indicates a greater loss of heat to outer space.The simulations Bio and CARB have the least negative net longwave heat flux, while HEAT has the highest negative heat flux, indicating that HEAT loses more heat than the other simulations.The higher heat loss in HEAT is due to the lowest SST and a reduced amount of greenhouse gases, specifically a low specific humidity (Fig. 6) and atmospheric CO 2 concentration (Fig. 4).The lower amount of greenhouse gases in the atmosphere permits a higher loss of heat outside the atmosphere.All the simulations considering phytoplankton light absorption, except CARB where the heat fluxes are prescribed, have a higher negative net longwave heat flux compared to Bio.This result is predictable because this biogeophysical mechanism is an additional heat source for the surface of the ocean, where air-sea heat exchanges occur.
The sensible heat flux depends on the atmospheric and oceanic temperature (Fanning and Weaver, 1996;Weaver et al., 2001).The sensible heat flux increases when the atmospheric temperature decreases and when the oceanic temperature increases.For the simulation HEAT, the sensible heat flux is the highest (Fig. 5c), because the atmospheric temper-ature is the lowest (Fig. 7).In contrast, the sensible heat flux is the lowest for the simulation BioLA, because the gradient of temperature between the ocean and the atmosphere is low.The sensible heat fluxes in HCorg and HCorgSI are close to the sensible heat flux of Bio because their air-sea temperature gradients are almost similar.
The global mean latent heat flux (Fig. 5d) depends mainly on the global mean precipitation (Weaver et al., 2001).The simulation Bio has the smaller latent heat flux due to the lowest precipitation in this simulation (Appendix G1).The latent heat fluxes in BioLA, HCorg, HCorgSI and HCorgSol are almost similar due to their almost similar precipitation values.The precipitation in HEAT is higher than in Bio, explaining the higher latent heat flux in HEAT.

Specific humidity and evaporation
The specific humidity and the evaporation in BioLA and HCorgSol are similar, and the same is true between the simulations HCorg and HCorgSI (Fig. 6).The specific humidity and evaporation are the lowest in Bio due to the lowest latent heat flux in this simulation.Including phytoplankton light absorption increases the latent heat flux and therefore increases the specific humidity and evaporation, which is consistent with Oschlies (2004) and Lengaigne et al. (2009).On a global scale, in BioLA the evaporation increases by 0.11 %, thus enhancing the specific humidity by 0.5 %.This latter value is lower than previous estimates where phytoplankton light absorption raises the specific humidity by 2 %-6 % (Patara et al., 2012).The different values between our estimates and Patara et al. (2012) come from the nonseasonal cycle in our model setup, changing the heat budget and therefore underestimating the specific humidity.Moreover, the specific humidity in HEAT is lower than in BioLA due to the lower latent heat flux in the simulation HEAT.The evaporation depends on several processes, and one of the most important is the humidity in the atmosphere, with https://doi.org/10.lower humidity leading to higher evaporation (Weaver et al., 2001).As a consequence, the evaporation is higher in HEAT than in BioLA.Furthermore, the specific humidity and the evaporation increase when the atmospheric temperature rises (Weaver et al., 2001).The specific humidity and evaporation are higher in CARB compared to BioLA because the surface atmospheric temperature is higher in CARB (Fig. 7).The specific humidity and evaporation in HCorg and HCorgSI are slightly lower than in BioLA because the latent heat flux in HCorg and HCorgSI is slightly lower.Once the CO 2 solubility factor is considered (simulation HCorgSol), the values of the specific humidity and evaporation are similar to the values in BioLA.This is not surprising because the heat fluxes between HCorgSol and BioLA are identical.

Surface atmospheric temperature
The difference in atmospheric properties between simulations lead to changes of the surface atmospheric temperature (Fig. 7).First of all, Bio has a low SAT because it does not include the additional heat source coming from the phytoplankton light absorption mechanism.The SAT in Bio is 9.31 • C, while the SAT in BioLA is 9.45 • C, which makes a global difference of 0.14 • C.This estimate is lower than previous estimates of 0.2-0.45• C (Patara et al., 2012;Asselot et al., 2021) due to the non-seasonal cycle in our model.
The lower SAT in HEAT compared to Bio is due to several reasons.Even if HEAT considers phytoplankton light absorption, we show that the SST in HEAT is lower than in Bio.For the SAT computation, the atmospheric CO 2 concentrations are identical between Bio and HEAT.Additionally, the specific humidity only increases by 0.27 % in HEAT compared to Bio.Therefore the greenhouse gas effect between these two simulations is rather similar.However, the global net longwave heat flux decreases by ∼ 0.2 W m −2 in HEAT due to the lower SST, leading to a cooling of the atmosphere.The combination of these different reasons explains the slightly lower SAT in HEAT compared to Bio.
For the simulation CARB, the concentration of greenhouse gases (atmospheric CO 2 and specific humidity) is higher than in Bio, while the air-sea heat fluxes are identical.As a consequence, more heat is trapped in the atmosphere, and the SAT increases by 0.71 • C compared to the reference run.
The sea-ice extent is different between HCorg and HCorgSI (Appendix F1), but the response of SAT is identical, indicating that with our model setup, the changes in sea-ice extent do not affect the SAT.The specific humidity and the atmospheric CO 2 concentration are slightly higher in HCorg and HCorgSI, leading to a small increase in SAT compared to Bio.
In HCorgSol the atmospheric CO 2 concentration and the specific humidity are higher than in Bio.However, the sen-  sible heat flux and the net longwave heat flux are lower in HCorgSol.Even if the greenhouse gas concentrations are higher, the reduced air-sea heat fluxes lead to a slight decrease in SAT in the simulation HCorgSol compared to Bio.

Limitations
Our study is designed to understand the climate pathways behind phytoplankton-induced atmospheric warming, but our model setup has limitations.Most notably, we do not consider a seasonal cycle in our study.Enabling seasonality would lead to larger seasonal increase in temperature, but it would also lead to larger seasonal decrease in CO 2 solubility.Therefore, we suppose that the heat pathway would not overrule the CO 2 pathway.Due to our new model setup with primary production allowed until the sixth layer of the model, the atmospheric CO 2 concentration of the simulations is lower than the pre-industrial level.Our quantitative estimates would be affected if the atmospheric CO 2 concentration of the reference run were higher, but we assume that the qualitative estimates would be very similar and the conclusions will not change.If the atmospheric CO 2 baseline is higher, the heat budget would also increase.As a consequence, the ocean would be warmer and the CO 2 solubility would decrease, increasing the importance of the CO 2 pathway in the phytoplankton-induced atmospheric warming.This study investigates the impact of short-lived seasonal organisms; thus having an annual mean approach underestimates the effect of phytoplankton light absorption on the climate system.The absence of a nitrogen cycle could have additional effects that are not included in our study.Phytoplankton light absorption warming low-oxygen regions causes an additional oxygen consumption, leading to hypoxia.As a consequence, we suggest that the denitrification would increase nitrogen fixation inducing a local increase in biomass.This increase in biomass would increase any pathway sensitivity of atmospheric CO 2 .Moreover, several studies considering a dynamic land component and focusing on phytoplankton light absorption report an increase in heat budget due to this biogeophysical mechanism (Anderson et al., 2007;Lengaigne et al., 2009).If a land model were to be included, we speculate that we would still find oceanic and atmospheric heating.However, the magnitude of changes might be smaller due to the uptake of CO 2 by vegetation, decreasing the atmospheric CO 2 concentration.Furthermore, the model does not include a temperaturedependency of iron bioavailability.According to previous experiments, a warming of the ocean decreases the bioavailability of iron (Liu and Millero, 2002).Phytoplankton light absorption increasing oceanic temperature might thus reduce the iron bioavailability.As a consequence, the limitation of phytoplankton growth by iron would increase, limiting the increase in chlorophyll biomass due to phytoplankton light absorption.Our study only considers two PFTs, and bring-ing in more PFTs would be an interesting complement to our findings.For instance, observations and modeling studies indicate that positively buoyant phytoplankton groups, such as cyanobacteria, are important to study phytoplankton light absorption (Sonntag and Hense, 2011;Paulsen et al., 2018;Wurl et al., 2018).Implementing these microorganisms to assess our research question could be a beneficial follow-up of our study.

Conclusions
For the first time, using the EcoGEnIE model (Ward et al., 2018), we compare the role of the air-sea heat and CO 2 fluxes and quantify their influence on the biologically induced atmospheric warming.We show that without any seasonality and with all the climate pathways included, the surface atmospheric temperature increases by 0.14 • C due to phytoplankton light absorption.As suggested by previous studies (Capone et al., 1998;Oschlies, 2004;Wetzel et al., 2006), phytoplankton light absorption changes the air-sea heat flux.Our results indicate that when only this air-sea interaction is considered, the atmosphere cools by 0.02 • C compared to a simulation without the biogeophysical mechanism.Moreover, when only the air-sea CO 2 exchange is considered, the atmospheric temperature increases by 0.71 • C. Clearly, our results indicate that the air-sea CO 2 exchange has a more important effect than the air-sea heat flux on the phytoplankton-induced warming of the atmosphere.With our model setup, the sea-ice extent and thickness slightly vary between simulations; therefore sea-ice processes hardly affect the air-sea CO 2 flux and thus the climate system.Moreover, including the solubility pathway changes the heat fluxes, specifically reducing the sensible heat flux and the net longwave heat flux compared to the reference simulation.As a consequence, this climate pathway has a negative effect on the atmospheric temperature.To conclude, phytoplankton light absorption influences the climate pathways at the ocean-atmosphere interface, particularly the air-sea CO 2 exchange that is important for the phytoplankton-induced atmospheric warming.For future climate studies, this work evidences that to capture the overall effect of climate-relevant mechanisms such as phytoplankton light absorption, the atmospheric CO 2 concentrations should evolve freely.
For future work, more studies with higher-complexity models are necessary to make quantitative assessments rather than qualitative assessments as in our study.Similar simulations must be conducted with a seasonal variation of the shortwave radiation to better understand the role of phytoplankton in the climate system.Moreover, a model with a dynamic atmosphere such as PLASIM-GEnIE (Holden et al., 2016) could be a good aspiration to complete our study.Indeed, previous studies evidence either an increased wind speed in subpolar regions (Patara et al., 2012) or enhanced atmospheric dynamics (Wetzel et al., 2006;Gnanadesikan Biogeosciences, 19, 223-239, 2022 https://doi.org/10.5194/bg-19-223-2022 and Anderson, 2009) due to phytoplankton light absorption.The increased wind speed with a dynamic atmospheric component could thus increase the air-sea CO 2 flux, reinforcing the importance of the CO 2 pathway in our study.Finally, implementing the new temperature-dependent remineralization scheme (Crichton et al., 2021;Armstrong McKay et al., 2021) would affect the biological pump and would be an extension to our findings.

Figure 1 .
Figure 1.Representation of the four different biologically induced pathways that affect the atmospheric temperature.(1) Marine biota via phytoplankton light absorption increases the SST, changing therefore the air-sea heat exchange and the atmospheric temperature.(2) Changes in SST also alter the solubility of CO 2 and its dissolved concentration.In turn, changes in dissolved CO 2 concentrations alter the air-sea CO 2 exchange and thus the greenhouse gas effect.(3) Phytoplankton light absorption modifies the marine biogeochemical cycles and particularly the export production of carbon.These changes in export production of carbon modify the dissolved CO 2 concentration and the greenhouse gas effect.(4) A warmer surface of the ocean can decrease the sea-ice extent.A reduction of sea-ice cover increases the air-sea CO 2 exchange area, changing the greenhouse gas concentrations.SAT: surface atmospheric temperature.SST: sea surface temperature.CaCO 3 : calcium carbonate.POC: particulate organic carbon.DOC: dissolved organic carbon.

Figure 2 .
Figure 2. Representation of the components of the EcoGEnIE model.The black arrows indicate the link between the different climatic components.

Figure 3 .
Figure 3. Sketch representing the climate pathways involved in the seven simulations (PLA: phytoplankton light absorption).Note that this figure is a simplification of Fig. 1; only the relevant pathways are represented.The name of the simulations are on the bottom left of each panel.The dashed arrows indicate the climate pathways prescribed.All the prescribed pathways are from the reference simulation Bio except the pathway between atmospheric CO 2 and SAT in CARB, which is prescribed from the simulation BioLA.

Figure 4 .
Figure 4. Atmospheric CO 2 concentrations (ppm) for the different simulations.The percentages represent the relative changes compared to the reference simulation Bio.

Figure 5 .
Figure 5. Global average of the different air-sea heat fluxes (W m −2 ) for the seven simulations.(a) Net shortwave radiation at the top of the atmosphere.(b) Net re-emitted longwave radiation.The net longwave radiation is negative, because it is dominated by the outgoing longwave radiation.(c) Sensible heat flux.(d) Latent heat flux.The color coding between the panels remains the same.The percentages represent the relative changes compared to Bio.

Figure 6 .
Figure 6.Evaporation (mm yr −1 ) and specific humidity (g kg −1 ) for the seven simulations.The percentages represent the relative changes compared to Bio.The name of the simulations are on the bottom left of each panel.

Figure 7 .
Figure 7. Sketch representing the surface atmospheric temperature (SAT) changes between the simulations and the reference run.On top is located the value of SAT change compared to Bio.The rest of the sketch is similar to Fig. 3.

Table 1 .
Sea surface temperature ( • C) and surface chlorophyll biomass (mgChl m −3 ).There is no value for the simulation CARB because we run the model with an uncoupled ocean-atmosphere setup.As a consequence, the surface PO 4 concentration is ∼ 0.216 µmol kg −1 in HCorg and ∼ 0.214 µmol kg −1 in HCorgSol.The higher surface PO 4 concentration leads to a larger surface chlorophyll biomass and higher SST in HCorg compared to HCorgSol.