Marine phytoplankton stoichiometry mediates nonlinear interactions between nutrient supply , temperature , and atmospheric CO 2

Marine phytoplankton stoichiometry links nutrient supply to marine carbon export. Deviations of phytoplankton stoichiometry from Redfield proportions (106C : 1P) could therefore have a significant impact on carbon cycling, and understanding which environmental factors drive these deviations may reveal new mechanisms regulating the carbon cycle. To explore the links between environmental conditions, stoichiometry, and carbon cycling, we compared four different models of phytoplankton C :P: a fixed Redfield model, a model with C :P given as a function of surface phosphorus concentration (P), a model with C :P given as a function of temperature, and a new multi-environmental model that predicts C :P as a function of light, temperature, and P. These stoichiometric models were embedded into a fivebox ocean circulation model, which resolves the three major ocean biomes (high-latitude, subtropical gyres, and tropical upwelling regions). Contrary to the expectation of a monotonic relationship between surface nutrient drawdown and carbon export, we found that lateral nutrient transport from lower C :P tropical waters to high C :P subtropical waters could cause carbon export to decrease with increased tropical nutrient utilization. It has been hypothesized that a positive feedback between temperature and pCO2,atm will play an important role in anthropogenic climate change, with changes in the biological pump playing at most a secondary role. Here we show that environmentally driven shifts in stoichiometry make the biological pump more influential, and may reverse the expected positive relationship between temperature and pCO2,atm. In the temperature-only model, changes in tropical temperature have more impact on the 1 pCO2,atm (∼ 41 ppm) compared to subtropical temperature changes (∼ 4.5 ppm). Our multi-environmental model predicted a decline in pCO2,atm of ∼ 46 ppm when temperature spanned a change of 10 C. Thus, we find that variation in marine phytoplankton stoichiometry and its environmental controlling factors can lead to nonlinear controls on pCO2,atm, suggesting the need for further studies of ocean C :P and the impact on ocean carbon cycling.


Introduction
The discovery of large-scale deviations of phytoplankton stoichiometry from the Redfield 50 ratio in the past decade (Martiny et al., 2013a(Martiny et al., , 2013bWeber and Deutsch, 2010) has significant consequences for our understanding of the biological carbon pump and global carbon cycling (Galbraith and Martiny, 2015;Moreno and Martiny, 2018). Traditionally, the biological pump is thought to be controlled by a combination of the vertical nutrient flux and nutrient utilization efficiency (Sarmiento and Toggweiler, 1984). Evidence that 55 elemental stoichiometry is variable thus adds a new dimension to the biological pump, and may lead to higher than currently expected carbon export in subtropical regions (Emerson, 2014;Tanioka and Matsumoto, 2017;Teng et al., 2014). Global carbon export has been estimated to range between 5 and 12 Pg C/year (Boyd and Trull, 2007;Henson et al., 2011), but these projections have yet to incorporate the environmental controls on 60 C:Pexport. Variation in C:Pexport from Redfield proportions can be linked to environmental conditions. There are two leading environmental parameters thought to control C:Pexport; nutrients, predominantly phosphate concentrations, and temperature. Galbraith and Martiny used a simple three-box model to show that variable stoichiometry driven by phosphate availability could enhance the efficiency of the biological pump in the low-65 latitude ocean (Galbraith and Martiny, 2015). In contrast, Yvon-Durocher and co-workers (2015) used a meta-analysis of global temperature and stoichiometric ratios to propose that C:P increased 2.6-fold from 0° C to 30° C. Thus, it is unclear if differences in nutrient supply, temperature, or some combination of them, control the global variation in C:P of plankton and exported material. 70 There are two important ingredients missing from published studies that could alter the interactions among phytoplankton stoichiometry, carbon export, and atmospheric pCO2 (pCO2,atm). The first is the presence of two distinct low-latitude biomes, namely the equatorial upwelling regions and the macronutrient-depleted subtropical gyres. In direct observations and inverse model analyses, these two biome types appear to have unique 75 elemental compositions (DeVries and Deutsch, 2014;Martiny et al., 2013a;Teng et al., 2014). Thus, in order to properly represent global variations in surface plankton C:P and carbon export, it is essential to separately model macronutrient-limited subtropical gyres and equatorial upwelling zones.
The second missing ingredient is that environmental factors beyond nutrient 80 availability may impact the elemental composition of surface plankton and C:Pexport. Temperature, irradiance, and nutrient concentrations are all important environmental factors, which influence the physiology and stoichiometry of phytoplankton. However, surveys of phytoplankton C:P are insufficient to distinguish the separate effects of each factor on C:P due to strong environmental covariance. Cellular trait-based models use 85 detailed studies of phytoplankton physiology to determine how phytoplankton cells should allocate their resources as a function of environmental conditions, allowing us to model the interactive influence of temperature, nutrient concentrations, and irradiance on C:P ratios (Clark et al., 2011;Daines et al., 2014;Shuter, 1979;Talmy et al., 2014;Toseland et al., 2013). Numerous physiological mechanisms have been proposed to explain variation in 90 phytoplankton stoichiometry, including growth rate (Sterner and Elser, 2002), photoacclimation (Falkowski and LaRoche, 1991;Geider et al., 1996;Geider, 2004, 2005), nutrient-limitation responses (Garcia et al., 2016;Goldman et al., 1979;Rhee, 1978), and temperature acclimation (Rhee and Gotham, 1981;Toseland et al., a trait-based model has revealed that differences in ribosomal content and cell size between warm-water, oligotrophic environments and cold-water, eutrophic environments are important mechanisms driving stoichiometric variation in the ocean (Daines et al., 2014). Thus, linking biome-scale variations in environmental conditions with a detailed trait-based model of phytoplankton resource allocation and elemental composition may 100 enable us to more fully explore interactions among multiple ocean environmental factors, the biological pump, and pCO2,atm.
Here, we create a five ocean circulation box model, incorporating the three major ocean biomes, to study the feedback effects of variable stoichiometry on carbon export and pCO2,atm. We will explicitly address the following research questions: (1) How does 105 environmental variability influence marine phytoplankton cellular allocation strategies and in turn the elemental stoichiometric ratio? (2) What are the effects of changing environmental conditions on stoichiometric ratios, carbon export, and pCO2,atm?, and (3) What is the influence of the environmental conditions among the three major surface biomes on carbon export and pCO2,atm? 110 2 Methods

Box Model Design
To quantify the feedbacks between phytoplankton stoichiometry, carbon export, and pCO2,atm, we formulated a five-box ocean circulation model of the phosphorus and carbon 115 cycles in the ocean coupled to an atmospheric box. The foundation of our model is based on the models introduced in Ito and Follows (2003) and DeVries and Primeau (2009). Phosphorus is used to represent the role of nutrient availability in controlling stoichiometry and C export. We chose this over N because on long time-scales, P is commonly considered the ultimate limiting nutrient whereas N is only limiting productivity 120 and export on short time-scales (Tyrrell, 1999). On long time-scales, nitrogen fixation/denitrification will presumably adjust the N inventory. Our modeling is focused on long term steady-state outcomes and we would like to avoid issues associated with modeling the N cycle (like getting N-fixation and denitrification rates correct). Thus, we chose to use P as a representative for nutrient availability at long-term steady-state 125 biogeochemical equilibrium. The model includes three surface boxes, each corresponding to one of the major biomes: the tropical equatorial upwelling regions (labeled T), the subtropical gyres (labeled S), and the high-latitude regions (labeled H) (Figure 1). We define the oligotrophic subtropical gyre regions where the mean annual phosphate concentration is less than 0.3 μM (Teng et al., 2014), with the remainder of the surface 130 ocean assigned either to box T or box H based on latitude. We use these assignments to calculate the baseline physical properties of each region, including mean annual averaged irradiance and temperature. The subsurface ocean is divided into two regions: the thermocline waters that underlies the subtropical gyres and the equatorial upwelling regions (labeled M), and deep waters (labeled D) (DeVries and Primeau, 2009

145
To simulate the global transport of water between boxes, our model includes a thermohaline circulation (Tc) that upwells water from the deep ocean into the tropics, laterally transports water into the subtropics and high-latitudes, and downwells water from the high-latitude region to the deep ocean. Surface winds produce a shallow overturning circulation (Tw), that transports water from the thermocline to the tropics and 150 then laterally into the subtropics. These circulations create teleconnections of nutrient supply in the surface ocean boxes. A bidirectional mixing term that ventilates the deep box directly through the high-latitude surface box (fhd) represents deep water formation in the Northern Atlantic region and around Antarctica (Sarmiento and Toggweiler, 1984). The parameters Tc, Tw and fhd are considered adjustable parameters, which we calibrate using 155 phosphate data from WOA13 (Garcia et al., 2014). In order to simulate the movement of particles, we included export fluxes (Pt, Ps, and Ph) of organic phosphorus out of each surface box. The conservation equations of phosphorus are as follows: where P represents the concentration of phosphorus at a specific box, a represents 0.25 remineralization, b represents 0.75 remineralization, and V represents the volume of the specified box. Our box model simulates [P], alkalinity and various forms of C; total carbon in the surface boxes is partitioned into carbonate, bicarbonate, and pCO2. The global mean [P] is 175 prescribed according to the observed mean value (Garcia et al., 2014). The export of carbon is linked to phosphorus export using the C:Pexport ratio. To quantify the breakdown of carbon into these components, we model the solubility pump, using temperature and salinity to determine the partitioning of inorganic carbon among total carbon within a box. The global mean alkalinity is prescribed according to the observed mean ocean values but 180 is also subject to transport (Sarmiento and Toggweiler, 1984). Our box model simulates alkalinity and total inorganic carbon, which are conserved tracers from which the speciation of inorganic carbon in seawater can be calculated. Biome specific salinity and temperature are used to prescribe the solubility constants of CO2 in seawater and the bromine concentration, which is taken to be proportional to salinity. CO2 cycles through 185 the atmosphere via the air-sea gas exchange fluxes (fah, fas, fat). We used a uniform piston velocity of 5.5 x 10 -5 m s -1 to drive air-sea gas exchange (DeVries and Primeau, 2009;Follows et al., 2002). Quantifying the atmospheric concentration of carbon satisfies: where C represents the concentration of total carbon in a specific box , Sol is the solubility constant in a specified box, calculated from temperature (temp) and salinity (sal). We calibrated our model parameters (Tc, Tw, fhd) so that the macronutrients were at similar average values compared to World Ocean Atlas 2013 dataset for each location. 195 We tested the sensitivity of modeled pCO2,atm to the fluxes Tc, Tw, and fhd and found that with Tc = 20 Sv and Tw = 5 Sv (values that allowed the model to match [P] and alkalinity), pCO2,atm was sensitive to the value of fhd (Sarmiento and Toggweiler, 1984). Guided by values previously used in the literature, we set fhd to 45.6 Sv (Table 1) but we also present results for the nutrient-only stoichiometry model at two extreme values of fhd (18 and 108 200 Sv) (Figure 2). The functional dependence of pCO2,atm with changing subtropical and tropical [P] for each extreme value of fhd was quite similar, though the value of pCO2,atm for the high fhd simulation was approximately twice that of the low fhd simulation (Figure 2). We found that our value of 45.6 Sv provides a modern pCO2,atm value. Although the focus of this study is to determine the impact of low latitude biogeochemistry on pCO2,atm, we point 205 out that at Redfield stoichiometry, pCO2,atm increases by 100 ppm when fhd increased from its default value 45.6 Sv to 108 Sv. For certain values of the parameters, the model produced excessive nutrient trapping in the thermocline. In order to dampen the nutrient trapping, we tuned the remineralization depth. As such, 25% of the total export is respired in the thermocline with the remaining 75% exported into the deep ocean, leading to a 210 better match between the modeled and observed [P] in the thermocline box.  (Sarmiento and Toggweiler, 1984) 3-300 (Toggweiler, 1999)

Stoichiometric Models
To quantify and understand the feedbacks between carbon export and pCO2,atm, we embedded four stoichiometric models into our five-box ocean circulation model. Each model differs according to its complexity and how much environmental information they 225 utilize. These are a static Redfield model that assumes that C:Pexport is constant across environmental conditions, a nutrient-only model that uses surface [P] to predict C:Pexport (from Galbraith and Martiny, 2015), a temperature-only model that uses T to predict C:Pexport (modified from Yvon-Durocher et al., 2015, and a multi-environmental model that uses light, T, and [P] to predict C:Pexport.

2.2.1Static Redfield Model
Our control model uses a static Redfield stoichiometry. The Redfield ratio is based on an average value of organic carbon to phosphorus of 106:1.

Nutrient-Only Model
The nutrient-only stoichiometric model expresses phytoplankton C:P as a function of the ambient phosphate concentration: where the parameters = 6.9 10 −3 μM -1 and [ ] 0 = 6.0 10 −3 were obtained by regressing the reciprocal of C:P onto [P] (Galbraith and Martiny, 2015).

Temperature-Only Model
The temperature-only stoichiometric model expresses phytoplankton C:P as a function of temperature: where the parameters = 0.037°⁄ and = 5.5938 (Yvon-Durocher et al., 2015). The 245 temperature-only model was created to determine the temperature responses of logtransformed C:P ratios centered at 15°C.

Multi-Environmental Model
We created a multi-environmental model which predicts how cell size, biomass allocations 250 to biosynthesis and photosynthesis, and C:P ratios vary with temperature, light levels, temperature, and phosphorus concentrations. The multi-environmental factor model was derived from a non-dynamic physiological trait-based model. We used a theoretical cellular-allocation trait model based on phytoplankton physiological properties that divides the 'cell' into several functional pools which represent cellular investments in 255 biosynthesis, photosynthesis, and structure, and a storage pool, which represents variations in the level of P-rich molecules such as polyphosphates (full model equations can be found on GitHub: https://github.com/georgehagstrom/-bg-2017-367-/blob/master/CP.m). The functional pools are composed of biological macromolecules such as ribosomes, proteins, carbohydrates and lipids. The model predicts the size of each 260 pool as a function of light, T, and [P]. The size of each functional pool is modeled by using subcellular resource compartments, which connect the fitness of a hypothetical phytoplankton cell in a given environment to its cellular radius and the relative allocation of cellular material to photosynthetic proteins, ribosomes, and biosynthetic proteins. We assume that real phytoplankton populations have physiological behaviors that cluster 265 around the strategy that produces the fastest growth rate in each environment (Norberg et al., 2001), and use the stoichiometry of this optimal strategy to model the elemental composition of cellular material ( Figure 1).
Phytoplankton can accumulate large reserves of nutrients that are not immediately incorporated into the functional components of the cell (Diaz et al., 2016;Mino et al., 1998;270 Van Mooy and Devol, 2008; Mouginot et al., 2015). This storage capability varies among phytoplankton species, and depends on the particular nutrient under consideration: the cost for storing physiologically relevant quantities of nutrients is low for nutrients with low quotas such as phosphorus, in comparison to nitrogen and carbon. Thus, the phosphorus storage is assumed highly plastic in comparison to carbon storage (Moore et al., 2013). 275 Further, we assume that each cell dedicates a fixed fraction of its biomass to carbon reserves, and focus our modeling efforts on the variability of the stored phosphorus pool. To predict the size of the storage pool, we assume a linear relationship between stored phosphorus and ambient environmental phosphorus levels and used statistical modeling of an oceanic C:P dataset  to calculate the constant of proportionality.

280
The result is a relatively simple model that both qualitatively and quantitatively predicts the variation of C:P in phytoplankton.
Phytoplankton physiology is modeled through allocations of cell dry mass to three distinct pools: structure (S(r)), biosynthesis (E), and photosynthesis (L). Allocations satisfy: where the variables S, E, and L represent the specific allocations of cellular biomass. The specific allocation of biomass to the cell membrane is inversely proportional to the cell radius � � (Clark et al., 2011), which accounts for the changing relative volume of the cell membrane with radius. The structure pool includes the cell membrane plus wall 295 and other components ( ), which are not related to photosynthesis or biosynthesis and is given by: In an environment specified by T, [P], and light level (I), the growth rate of a cell using a given strategy is the minimum of the following growth rates: Here is determined by the specific rate of protein synthesis, is determined by the specific rate of carbon fixation, and is determined by the specific rate of phosphorus uptake, or: We assume that part of the energy captured by a cell via photosynthesis is used for 305 maintenance ( ), whereas the rest is used to drive the synthesis of new macromolecules ( ), so that a cell growing at rate is in energy balance. The efficiency of biosynthesis and the carbon cost of maintenance are functions of T, whose dependence is modeled using 10 = 2.0 (Van Bogelen and Neidhardt, 1990;Broeze et al., 1978;Shuter, 1979). Uptake is regulated by a Monod function with kinetic parameters depending on the radius 310 through the allometric scaling relationships derived from measurements of phytoplankton populations (Edwards et al., 2012): This use of allometric scaling relationships departs from the conventions adopted by Shuter (1979) or Daines et al. (2014), who assume that uptake rates are diffusion-limited.

315
The phosphorus quota for functional elements of the cell (thus not including any storage) is determined by the allocation to biosynthesis and the percentage of cellular dry mass allocated to DNA: Here, we assume that there is no contribution to the functional-apparatus P quota from 320 phospholipids, which instead are merged with storage molecules. This differs from Daines et al. (2014), who assumes that phospholipids occupy 10% of the cell by mass. Phytoplankton can substitute sulfoquinovosdiaglycerol (SQDG) for phospholipids in their cell membranes under low P conditions (Van Mooy et al., 2009). Similarly, P storage molecules are also regulated by P availability. Thus, we treat phospholipids and P-storage 325 as one pool. The function is the cellular response to light levels, and is chosen to capture the effects of both electron transport and carbon fixation on photosynthesis, and is closely related to a previous model (Talmy et al., 2013). This prior model included four compartments: electron transport, carbon fixation, photoprotection, and biosynthesis. It 330 was found that photoprotection allocation was not a large or greatly changing component of their allocations. We therefore do not include this within our model due to its high complexity with little qualitative results. Our biosynthesis was also separately parametrized.
The decomposition of photosynthesis into light harvesting and carbon fixation 335 components is critical, and makes our model predictions agree much better with experiments studying the variations of C:P or N:P ratios with irradiance. Models that do not have this decomposition predict too large a decrease in cellular allocations to photosynthesis at high-light levels. In a two-compartment model, increases in allocations to carbon fixation cause the overall allocation to light harvesting to have a more mild 340 decrease. The two-compartment treatment also seems more physiologically realistic than a 1-compartment treatment, which only models photosynthetic pigments. Thus, we used the functional forms and parameters that were derived (experimentally) previously for carbon fixation and light harvesting (Talmy et al., 2013). Our model interprets light harvesting allocation, , as being composed of proteins 345 dedicated to carbon fixation ( 1 ), such as RuBisCO, and proteins dedicated to light harvesting ( 2 ), such as photosynthetic pigments. The rate of photosynthetic carbon fixation is a function of the allocations to each of these, which satisfy 1 + 2 = . The relative allocations together determine the overall photosynthetic rate: 350 For a given I and , there is a pair of values � 1,opt , 2,opt � that maximize the photosynthetic rate . We estimate the photosynthetic rate ( , ) under the assumption that cells assume the optimal allocations to carbon fixation and electron transport. Our model departs from the models developed by Shuter (1979) and Daines et al. (2014), which assume that energy acquisition is a linear function of light levels leading to functional 355 responses linearly proportional to the cellular investment in light harvesting proteins. We model photosynthesis as having a Q10=1, which is consistent with physiological studies going back to Shuter (1979) that suggest that photosynthetic efficiency does not depend on temperature over physiologically relevant ranges. The discrepancy between photosynthetic and biosynthetic temperature dependence has traditionally been explained 360 by referring to the differences in the chemistry and physics of the two processes. The electron transport chain relies on quantum mechanical processes, which are unaffected by variations in temperature in a physiologically relevant range (Devault, 1980). Required maintenance respiration rates are also modeled as having a Q10=2.0 (Devault, 1980). Required maintenance respiration rates are also modeled as having a Q10=2.0. We model 365 the phytoplankton community residing in a given environment by assuming it consists solely of the phytoplankton type using the highest growth rate strategy in that environment. This strategy is found by solving for the values of and and that make = = = (16) We will now show that under two assumptions that will be true in nearly any realistic 375 situation, a strategy maximizing always exists, is unique, and satisfies = = = (Figure 4). The function is a function of the chosen strategy ( , ), and it is an increasing function of and decreasing function of . The first assumption is that light levels are sufficiently high that there exists some such that ( , 0) > 0, which means that light is sufficient for some phytoplankton to be able to overcome maintenance costs. The 380 function is a monotonically decreasing function of both and . As there is a non-zero amount of contained in the structure pool, and because uptake rates decline to zero with , there will be some at which ( , 0) > 0. The second assumption is that < , which will be true for most realistic values of the light level. We note that for fixed , is a monotonically decreasing function of . Since none of , , or have 385 critical points, the function can only have a maximum at places where two or more of , , and are equal, or at the boundaries of the strategy space. On the boundaries of strategy space, = 0 or = 0 so that ≤ 0. We can exclude the boundary and focus on places where two or more of , , and are equal. We define two curves, one on which = , and the other on which = . The curve for which = begins at the point 390 = and can be described by a monotonically increasing function = ( ) on the interval [ , ∞]. This curve exists because = 0 when = 0, > 0 when = 0 and < , and < 0 when = 1 − ( ) − = 0, so that there is always a solution to = for fixed > . To see that the curve is an increasing function of , consider the function ( , ) = − and apply the chain rule to the equation ( ( ), ) = 0 to find 395 that along the curve E=g(r) : We consider the terms in equation 17 carefully. The function is an increasing function of because is independent of and because is an increasing function of (for a fixed 400 investment in biosynthesis, a larger radius leads to a greater investment in photosynthesis and greater photosynthetic growth rate). Thus, the numerator of equation 17 is negative. The function is a decreasing function of because is a decreasing function of (greater investments in biosynthesis at fixed radius lead to smaller investments in photosynthesis) and is an increasing function of . Thus the denominator of equation 17   405 is negative, and the quotient on the right hand side is positive, so ′( )is positive and describes an increasing curve. By similar logic, we can define a curve ℎ( )that solves the equation ( . Therefore, the only place where can have a maximum is at the place where ( ) and ℎ( ) intersect. This is the strategy that leads to equality of all the growth rates. We refer to this strategy, as a function of environmental conditions, as � ( , , ), ( , , ), ( , , )�. Using this strategy, we can predict the stoichiometry of the functional components of the phytoplankton population in a given environment. 420 We assume that real phytoplankton populations cluster near the optimal strategy in the local environment (Norberg et al., 2001): ( , ) = argmax ( , ) .
(18) For all values of environmental parameters used in this study, the unique maximum of the growth rate occurs for the set of parameter values that lead to co-limitation by nutrients, photosynthesis, and biosynthesis, analogously to the predictions of Klausmeier and co-425 workers (2004). The optimal strategy determines the model prediction of the C:P of functional components in a given environment by taking the quotient of the carbon and phosphorus quotas.  Shuter, 1979) apparatus at T 0 =25 Q 10,E Q 10 of biosynthetic apparatus 2.0 (Shuter, 1979) Φ M0 Specific carbon cost of maintenance at T 0 =25 10 −3 hr −1 (Shuter, 1979) Q 10,M Q 10 of maintenance 2.0 - (Shuter, 1979) Q10,P Q10 of photosynthesis 1.0 (Shuter, 1979) ΦS Carbon cost of synthesis 0.67 - (Shuter, 1979) (Sterner and Elser, 2002) The carbon quota is calculated as: Here we see the contributions of carbon contained in both functional and storage pools, the latter of which are assumed to occupy a fixed fraction of the cell independent of the environment (but linked to cell size). Measurements of cellular P partitioning indicate that the ribosomal RNA can 435 sometimes contribute only 33% of the total P quota (Garcia et al., 2016). The additional phosphorus includes membrane phospholipids and storage compounds, each of which can be up-or down-regulated in response to phosphorus availability in the environment. To model this phenomenon, we assume the existence of an additional stored P pool, whose size is a linear function of environmental P, or: where is determined by the best fit to the Martiny et al. (2014) data. Our model then predicts C:P as: The model parameter is calculated by minimizing the residuals of the P:C ratio predicted by Eq.19 in comparison to the global data-set on particulate organic matter stoichiometry 445 compiled by Martiny and others (2014). To maintain consistency with the linear regression model of Galbraith and Martiny (2015), we restrict the dataset to observations from the upper 30 meters of the water column containing particulate organic phosphorus and carbon concentrations of greater than 5 . Observations from the same station and the same day, but at different depths in the water column are averaged together. The P:C ratio 450 of the functional apparatus is calculated using irradiance, T, and [P] data from the World Ocean Atlas (Garcia et al., 2014;Locarnini et al., 2013; oceancolor.gsfc.nasa.gov/data/10.5067/AQUA/MODIS/L3B/PAR/2014/), which are used to estimate environmental conditions at the location and date of particulate organic matter measurements. Light levels are computed by averaging irradiance over the top 50 meters 455 of the water column, assuming an e-folding depth of 20 meters. Linear regression determines = 2500 M −1 which fits the data with an 2 = 0.28. All parameters for the model are listed in Table 2.

460
To address how changing environmental conditions affected stoichiometric ratios, carbon export, and pCO2,atm we performed two tests; a change in nutrients and a change in sea surface temperature. These tests allowed us to observe how the relationships between environmental conditions, carbon export and pCO2,atm, depend on the mechanisms responsible for stoichiometric variation in the ocean. In order to account for the effects of particulate inorganic carbon (PIC) export, we multiply model predicted C:Pexport by 1.2, consistent with previous studies (Broecker, 1982;Sarmiento and Toggweiler, 1984).
The first set of numerical experiments examined the sensitivity of pCO2,atm to nutrient availability in the tropical and subtropical boxes for each of the three stoichiometric models. We varied tropical [P] from 0.15 to 1.5 μM and subtropical [P] from 470 1x10e -3 to 0.5 μM by adjusting the implied biological export and determined the equilibrium pCO2,atm values.
The second set of experimental tests was done to quantify how temperature modifies carbon export and pCO2,atm for each stoichiometric model. Temperature influences carbon cycling in two ways within our model: through the solubility of inorganic carbon in 475 seawater and through changes in phytoplankton stoichiometry within the temperatureonly and multi-environmental models. Due to the well-known effects of temperature on CO2 solubility, it is generally predicted that there is to be a positive feedback between pCO2,atm and temperature mediated by declining CO2 solubility at high temperatures. To study the relative strengths of the temperature solubility feedback and the temperature 480 regulation of C:P feedback, we performed a numerical experiment in which we varied the sea surface temperature by five degrees in either direction of modern sea surface temperature. This represents a plausible range of variation under both ice-age and anthropogenic climate change scenarios. We varied tropical temperature from 21° to 31°C and subtropical temperature from 19° to 29°C, determining equilibrium pCO2,atm values for 485 combinations of temperature conditions.

Results
To quantify the linkages between phytoplankton physiology, elemental stoichiometry, and ocean carbon cycling, we divide our results into two parts. The first is a direct study of the 490 stoichiometric models, comparing their predictions about the relationship between stoichiometry and environmental conditions, and in the case of the trait-based model, illustrating how cellular physiology is predicted to vary across these conditions. In the second part, we show how variable stoichiometry influences carbon export and pCO2,atm, under changing phosphorus concentrations and temperature. Within these results, we 495 distinguish the influence or lack thereof off the three distinct biomes; in particular the equatorial upwelling regions and the macronutrient depleted subtropical gyres.

Multi-environmental and physiological controls on plankton stoichiometry
Our multi-environmental model captured several major mechanisms hypothesized to be 500 environmental drivers of C:P ratios including a temperature dependence of many cellular processes, a link between growth rate and ribosome abundance, and storage drawdown during nutrient limitation. The predicted relationship between environmental conditions and C:P can be understood through the environmental regulation of three factors: (i) the balance between photosynthetic proteins and ribosomes, (ii) the cell radius and associated 505 allocation to structural material, and (iii) the degree of phosphorus storage. Our model predicted that for an optimal strategy, specific protein synthesis rates will match specific rates of carbon fixation. Thus, the ratio of photosynthetic machinery to biosynthetic machinery is primarily controlled by irradiance and temperature. Increases in light levels lead to higher photosynthetic efficiency, higher ribosome content, smaller cells (due to a 510 lower requirement for photosynthetic machinery), and lower C:P ratios ( Figure 5). The response of C:P to light levels predicted by our model was muted in comparison to other subcellular compartment models because we separately modeled electron transport and carbon fixation (Talmy et al., 2013), and our predictions were consistent with the weak relationship between irradiance and C:P (Thrane et al., 2016) (Figure 5A).
Increases in temperature increase the efficiency of biosynthesis, but not photosynthesis. Therefore elevated temperature lead to a reduced ribosome content relative to photosynthetic proteins and higher C:P ratios ( Figure 6A). There leads to a nonmonotonic, concave relationship between temperature and cell size, which is due to a subtle interaction between biosynthesis efficiency (which varies greatly with temperature), 520 maintenance costs, and size dependent uptake rates.
Nutrient concentrations do not affect the ratio of biosynthetic to photosynthetic machinery but positively relate to both P storage and cell radius. Cell radius directly influences the specific rate of nutrient uptake, and indirectly biosynthesis and photosynthesis as the cell membrane and wall affects the space available for other 525 investments. This effect is pronounced in oligotrophic conditions ([P] < 100nM). Here, cell radius declines below 1μm resulting in decreasing allocations to both photosynthesis and biosynthesis and elevated C:P ratios. Higher values of the cell radius are observed in nutrient rich conditions. P concentrations also influenced C:P through the direct control of P storage. We 530 plotted the relative contribution of the storage compartment and the functional compartment to the P quota, as a function of environmental conditions. The impact of the residual pool on the overall size of the P pool is heavily dependent on environmental conditions, varying from a minimum of close to 0% to a maximum of just under 50%. In the vast majority of the parameter range considered here, the contribution of the residual pool 535 is much more modest (10-20%). High values occur when phosphorus is available and the temperature is high. In these conditions, ribosomal contributions are decreased but the residual contribution is high. In cold water, high P ecosystems, the residual contribution is approximately 25%, and in oligotrophic ecosystems it is close to 0. Thus, C:P was predicted to be a decreasing function of [P] with two distinct regimes: a moderate sensitivity regime 540 for [P] above 100nM, and a high sensitivity regime for [P] below 100nM.

Figure 5: Influence of phosphate concentration and irradiance on cellular stoichiometry and cellular traits, at a constant T = 25 °C. A) Cell radius (r). B) P storage allocation. C) Biosynthesis allocation. D)
545 Photosynthesis (L) allocation. E) The C:P ratio. As irradiance increases, there is a tendency towards greater allocation to biosynthesis and lesser allocation to photosynthesis, which leads to lower C:P ratios. When phosphorus is very low, there is a large decrease in both biosynthesis and photosynthesis allocations due to the large relative allocation to the cell membrane. C:P ratios are inversely proportional to phosphorus concentration, driven by an increase in luxury storage and ribosomal content as P increases. We next used the outcome of the trait model as a multi-environmental model to predict C:P ratios in the modern ocean based on annual mean light, T, and [P]. Our predictions reproduced the global pattern (Martiny et al., 2013b) with C:P ratios above the 560 Redfield ratio in subtropical gyres and C:P ratios below the Redfield ratio in equatorial and coastal upwelling regions and subpolar gyres ( Figure 7A). Additionally, our model also reproduced basin-scale stoichiometric gradients among similar biomes in each ocean, predicting the highest C:P ratios in the western Mediterranean Sea and the western North Atlantic Subtropical Gyre, and somewhat elevated C:P ratios in the South Atlantic 565 Subtropical Gyre as well as the North and South Pacific Subtropical Gyres.

Figure 7: Predicted C:P ratios in the global ocean in differing climatic regimes. A) C:P ratio under modern
ocean conditions. Large differences in C:P are predicted between distinct types of ocean biome, with low C:P in 570 equatorial upwelling regions and subpolar gyres, and high C:P in subtropical gyres. Regional differences between biomes of similar type are observed as well, with the low phosphorus Atlantic having a higher C:P than the Pacific. B.) C:P ratio under cooling temperature conditions (-5°C from the modern ocean). C) C:P ratio under warming temperature conditions (+5°C from the modern ocean). Each 5 degree change leads to a shift of 15% in the mean C:P ratio of organic matter.

575
To study the potential impact of sea surface temperature on phytoplankton resource allocation and stoichiometry, we used our multi-environmental model to predict C:P in ocean conditions both five degrees colder (Cooling environments) and warmer (Warming environments) than the modern ocean. According to our model, a five-degree increase (or 580 decrease) in sea surface temperature would cause a 15% rise (or fall) in C:P ratios ( Figure  7). This sensitivity suggested that the relative effect of T on biochemical processes could have large implications for biogeochemical cycles, making it important to determine the relative importance of physiological mechanisms regulating C:P ratios. We compared the multi-environmental model to the predictions made by two other models: the nutrient-only model used by the Galbraith and Martiny model (2015), and our temperature-only model modified from Yvon-Durocher and co-workers (2015). These two 595 models also successfully predicted the qualitative pattern of stoichiometric variation in the ocean, but were unable to replicate the full range of variation observed in the data ( Figure  8). In particular, there were mismatches in the North Atlantic Subtropical Gyre and the Southern Ocean, where the C:P ratio is at the extreme ( Figure 8A, B). The nutrient-only model had a tendency to predict lower C:P ratios than the multi-environmental model in 600 warm tropical and subtropical waters, and predict higher C:P ratios in cold waters ( Figure  8A). This difference is driven by the T sensitivity of biosynthesis in the multi-environmental model, leading to increasing C:P in all warm water regions and decreasing C:P in cold water regions ( Figure 8C ). The multi-environmental model predicted a wider range of C:P in the ocean. The temperature-only model overall had higher C:P ratios 605 globally compared with the multi-environmental model ( Figure 8B) but suggested lower C:P in the gyres and higher C:P in high latitudes, especially in the Southern Ocean ( Figure  8D). 610 We next quantified the impact of nutrient availability in the tropics and subtropics on stoichiometry, carbon export, and pCO2,atm ( Figure 9A-L). Using a constant Redfield model (or the temperature-only model), we replicated the previously observed approximately linear relationship between surface [P] and pCO2,atm (equivalent to how pre-formed [P] will influence pCO2,atm) (Ito and Follows, 2003;Sigman and Boyle, 2000). We found that [P] 615 drawdown in the subtropical box had a greater impact on carbon export, since export from the high-latitude box was not enhanced by the [P] supply from the subtropical box ( Figure  9A, D, G). In the Redfield model, pCO2,atm appeared to be much more sensitive to subtropical [P] than tropical [P], which was partially due to enhanced carbon export in the subtropical box and partially due to the larger surface area of the subtropical box (implying 620 a greater potential for CO2 exchange) ( Figure 9J).

Impact of nutrient availability on carbon export and atmospheric pCO2
In contrast to the predictions made using Redfield stoichiometry, when we used the nutrient-only model for phytoplankton stoichiometry, we observed a non-linear relationship between surface [P] and pCO2,atm ( Figure 9B, E, H, K). At fixed tropical [P], there was a strong relationship between subtropical [P] drawdown, export, and pCO2,atm in 625 accordance with the findings of Galbraith and Martiny (2015) (Figure 9B, E,H). The total decline in pCO2,atm as subtropical [P] declined from 0.4 μM to 1x10 -3 μM could be more than 60 ppm, which was more than twice the decline that occurred in the fixed stoichiometry experiment ( Figure 9K). We found a non-linear monotonic relationship between tropical [P] and pCO2,atm: when tropical [P] was high, declines in tropical [P] led to 630 lower carbon export and increased pCO2,atm. However, this trend reversed when tropical [P] was further drawn down ( Figure 9K). The counter intuitive decline in pCO2,atm with higher export from tropics was driven by a teleconnection in nutrient delivery between the subtropical and tropical boxes. Increases in export in the tropical box due to [P] drawdown decreased the supply of [P] to the subtropics, which led to a decrease in the more efficient 635 (higher C:P) subtropical export. Thus, the nutrient-only model predicted a greater decrease in subtropical export than the increase in tropical export.
The multi-environmental model also predicted a non-linear relationship between P draw down, carbon export, and pCO2,atm. However, the pattern was somewhat distinct from that of the nutrient-only model results ( Figure 9C, F, I, L). First, subtropical [P] drawdown 640 had a nonlinear relationship with pCO2,atm: when subtropical [P] was high, declines in tropical [P] led to slight declines in pCO2,atm, and when subtropical [P] is low, small declines in tropical [P] lead to large declines in pCO2,atm. This intensification of the relationship between subtropical [P] and pCO2,atm was due to the nonlinear relationship between [P] and C:P predicted by the trait-based model ( Figure 9I). The multi-environmental model 645 predicted extremely high tropical export, but only when [P] was lower than 0.05 μM ( Figure 9C, F, I). Second, the effect of tropical [P] levels on pCO2,atm was strongly modulated by subtropical [P], reversing from a negative to a positive relationship as subtropical [P] declines ( Figure 9I, L). The difference between the nutrient-only model and the multienvironmental model arose because the multi-environmental model incorporated a 650 temperature impact on resource allocation and elemental ratios. Although we were not varying temperature in these experiments, we did represent regional temperatures differences between the different boxes. The result is that a large stoichiometric contrast between the tropical and sub-tropical regions only arose when there was a large difference in nutrient levels between the two regions (Fig. 9L). However, both the nutrient-only model 655 and the multi-environmental model predicted that carbon export and pCO2,atm were sensitive to the interaction between regional nutrient availability and C:Pexport.

660
Columns correspond to type of stoichiometry; Redfield (Left), nutrient-only (Middle), and multi-environmental model (Right). Rows correspond to either tropical carbon export (A through C), subtropical carbon export (D through F), total carbon export (G through I) or atmospheric pCO2 (J through L). The grey points represent where pCO2,atm was calculated, between spaces are interpolated.

Interactive effect of temperature on stoichiometry, carbon export and atmospheric pCO2
We next quantified the impact of sea surface temperature (SST) in the tropics and subtropics on C:Pexport, carbon export, and pCO2,atm ( Figure 10A-D). The Redfield model predicts that increases in temperature lead to a decline in the solubility of CO2 in seawater 670 and consequently an increase in pCO2,atm from 288 to 300 ppm (Δ pCO2,atm = 12) ( Figure  10A). This feedback was present with the same strength in the nutrient-only model (with no T dependence on C:P), in which pCO2,atm ranged from 268 to 280 ppm (Δ pCO2,atm = 12) ( Figure 10B).
In contrast to the Redfield and nutrient-only models, the temperature-only model predicted a negative linear relationship between pCO2,atm and tropical sea surface T and a positive linear relationship between pCO2,atm and subtropical sea surface T ( Figure 10C). The decline in pCO2,atm with tropical SST was driven by an enhancement of export due to increased C:P at higher temperatures ( Figure 11). At 5°C below modern ocean temperature, the model predicted C:P in the tropics was 131 and subtropical was 121, resulting in a 680 pCO2,atm of 305 ppm. At 5°C above modern ocean temperature, the model predicts a C:P ratio in the tropics of 189 and C:P ratio of 175 in the subtropics, resulting in a pCO2,atm of 263 ppm. Tropical SST had more impact with Δ pCO2,atm = 41 ppm compared to subtropical SST with a ΔpCO2,atm ranging from 4 to 5 ppm ( Figure 11). Similar to the temperature-only model, the multi-environmental model predicted a 685 negative linear relationship between pCO2,atm and tropical SST and a positive linear relationship between pCO2,atm and subtropical SST ( Figure 10D). The decline in pCO2,atm with tropical SST was driven by an enhancement of export due to increased C:P at higher Ts ( Figure 11). In the subtropical region, the expected increase in export was mitigated by a decline in solubility. At 5˚C below modern ocean temperature, the trait-based model 690 predicted that C:P in the tropics was 147 and that C:P in the subtropics was 155, resulting in an increase of pCO2,atm to 279 ppm ( Figure 11). Variation in tropical SST over a 10˚C span led to a significant decline in pCO2,atm, with a Δ pCO2,atm of approximately 46, and tropical C:P ranging from 147 to 210 ( Figure 11). Because the subtropical box has a large surface area, the decrease in surface CO2 solubility at high temperatures is sufficient to overcome the increase in export due to higher C:P leading to a positive relationship between pCO2,atm and subtropical temperatures.

710
Here, we found that variable stoichiometry of exported organic material moderates the interaction between low-latitude nutrient fluxes and ocean carbon cycling. A full connecting circulation allows for complete movement of nutrients between ocean regions resulting in strong linkages between nutrient supply ratios and cellular stoichiometric ratios (Deutsch and Weber, 2012). It has been shown that the inclusion of an oceanic 715 circulation connecting high and low-latitude regions results in a feedback effect between high-latitude nutrient export and relative nutrient supply in low-latitudes (Sarmiento et al., 2004;Weber and Deutsch, 2010). Together, the inclusion of lateral transport between ocean regions and of deviations from Redfield stoichiometry within our model led us to predict the existence of strong teleconnections between the tropics and the macronutrient 720 limited subtropics. The degree of nutrient drawdown in the tropics had a strongly nonmonotonic relationship with pCO2,atm because this drawdown influenced both nutrient supply to the subtropics and tropical C:P. The idea of biogeochemical teleconnections has been proposed before, but we found that variations in stoichiometry greatly enhance the importance and strength of such linkages (Sarmiento and Toggweiler, 1984). Thus biomescale variations in phytoplankton elemental stoichiometry may change the sensitivity of the carbon pump to other phenomena that regulate patterns of nutrient drawdown. We also see that the degree of nutrient drawdown had a strong impact on predicted (and observed) C:P leading to highly non-linear controls on pCO2,atm, whereby increased export in the tropics counter intuitively leads to increasing pCO2,atm. Large-scale gradients in 730 stoichiometry can alter the regional efficiency of the biological pump: [P] supplied to high C:P regions leads to a larger export of carbon than [P] supplied to low C:P regions. This lends an important role to details of ocean circulation and other processes that alter nutrient supply and phytoplankton physiological responses in different surface ocean regions. Therefore, biome-scale variations in phytoplankton elemental stoichiometry can 735 lead to a fundamental change in the partitioning of carbon between the atmosphere and the ocean.
We have created a box model to simulate the impact of the low latitude stoichiometric ratios, its environmental controlling factors, and the relationships to pCO2,atm. Low latitude phosphorus concentrations can be set in one of two fashions; 740 through iron limitation and through nutrient supply. Here we will briefly discussion how iron limitation would play a significant role on phosphorus concentrations and associated C:P. The biogeochemical functioning of tropical regions are commonly influenced by iron availability in such a way that macronutrients cannot be fully drawn down by phytoplankton (Coale et al., 1996;Moore, 2004;Raven et al., 1999). The degree of nutrient drawdown has a strong impact on predicted (and observed) C:P. This environmental control on C:P could lead to highly non-linear controls on pCO2,atm whereby increased iron availability lead to increased [P] draw down and export in the tropics. However, as we have shown this may lead to increasing rather than commonly assumed decreasing pCO2,atm. This link between iron and export would differ in the subtropics, where iron is thought to 750 stimulate nitrogen levels through nitrogen fixation. This would result in elevated phosphate draw down, higher C:P and higher export. Thus, iron availability may play a complex role depending on whether there is an increased delivery in upwelling zones (leading to a potential declining global C export) or in the subtropical gyres (leading to a potential increase in global C export).

755
Past studies using box models have found pCO2,atm to be insensitive to low-latitude nutrients (Follows et al., 2002;Ito and Follows, 2003;Sarmiento and Toggweiler, 1984;Toggweiler, 1999). This phenomena was explored by DeVries and Primeau (2009), who showed that the strength of the thermohaline circulation is the strongest control on pCO2,atm, and that changes in low-latitude export has a minor impact. Unlike our study, such 760 earlier work relied on a uniform Redfield stoichiometry. However, we find that when stoichiometric variation is included, carbon export and pCO2,atm become dependent on details of low-latitude processes.
It is important to recognize that a five-box model is an incomplete description of ocean circulation, and is here used to illustrate important mechanisms, not to make precise 765 quantitative predictions. In order for our model to adequately reflect important features of the carbon and phosphorus nutrient distributions, we had to carefully select the values of the thermohaline and wind-driven upper ocean circulations that lead to reasonable nutrient fluxes and standing stocks. The value of thermocline circulation, Tc, has been calibrated in different box models to range from 12 to 30 Sv (DeVries and Primeau, 2009;770 Galbraith and Martiny, 2015;Sarmiento and Toggweiler, 1984;Toggweiler, 1999). Variations in the thermohaline circulation influence the abundance of nutrients in different boxes. Depending on the strength of this circulation, our model accumulated nutrients in the thermocline box and we tuned this parameter to most accurately mimic nutrient variation across ocean regions. Other caveats relates to our choice of the wind driven 775 overturning circulation and the two-way flux values. Similar to the circulation values, a wide range of two-way flux values have been used in the literature. We therefore performed sensitivity experiments to find the best value for our full model set-up but the qualitative trends observed are insensitive to the choice of such fluxes.
Nutrient availability and temperature have been alternatively proposed as drivers 780 of variation in stoichiometric ratios in the global ocean, and the strong statistical correlation between temperature and nutrients throughout the ocean has prevented identification of the relative importance of each factor (Martiny et al., 2013b;Moreno and Martiny, 2018). We see that although temperature regulation of C:Pexport can influence pCO2,atm, this regulation is strongly dependent on the detailed physiological control 785 mechanism and also generally diverge from expectations based on the solubility pump. The decrease in surface CO2 solubility at elevated temperature is sufficient to overcome the increase in export due to higher C:P leading to a positive relationship between pCO2,atm and subtropical temperatures. It is important to point out that the relative importance of the two competing effect depends critically on the physical circulation of the ocean. Predicted 790 increases in stratification are often invoked as a mechanism that would decrease the vertical supply of nutrients, which one might think would further compensate for the effect of higher C:P. However, the strength of the biological pump in the subtropics is also influenced by lateral transport of nutrients (Letscher et al., 2015) as such we argue that it is unclear if you should expect increasing, unchanged, or decreasing C export in low latitude 795 regions with ocean warming and stratification. Similarly, it is unclear how increases in stratification might affect the strength of the solubility pump. The sensitivity of pCO2,atm to changes in subtropical surface temperatures depends critically on the volume of the ocean ventilated from the subtropics, i.e. on the volume of the thermocline box in our model. How this volume might change in response to a warming world is a complicated dynamical 800 problem that is beyond the scope of the present work. Our results do not identify whether temperature or nutrient concentrations is the most important driver of phytoplankton C:P, but do suggest that the physiological effect of temperature could be important for ocean carbon cycling. Both the temperature-only and multi-environmental models predict that temperature increases enhance tropical export, 805 causing substantial decreases in pCO2,atm with temperature. This relationship is the reverse of that predicted by the nutrient-only and Redfield models, and represents a sizable potential negative feedback on carbon cycling. The multi-environmental model also predicted that C:P responds in a nonlinear fashion to [P], with significantly increased sensitivity in highly oligotrophic conditions. Thus, a deeper understanding of the 810 physiological mechanisms regulating phytoplankton C:P ratios are key to understanding the carbon cycle.
Our derivation of the multi-environmental model relies on several important assumptions. The growth rate in the multi-environmental model is determined by a set of environmental conditions and quantified by the specific rate of protein synthesis, carbon 815 fixation, and phosphorus uptake. The effect of growth rate on stoichiometry will likely be dependent on whether light, a specific nutrient, or temperature controls growth (Moreno and Martiny, 2018). The value of specific values of Q10 leads to uncertainty in our multienvironmental model because of the range of possible values is highly dependent on the cell or organism being tested. In a study examining Q10 of various processes within the cell, 820 it was found that the Q10 of photochemical processes ranged from 1.0 to 2.08, and for carboxylase activity of RuBisCO to be 2.66 (Raven and Geider, 1988). In addition to the high uncertainty between Q10 values, there is high ambiguity associated with cellular inorganic P stores (e.g., polyphosphates and phospholipids) (Kornberg et al., 1999). P storage, such as polyphosphates, can serve as both energy and nutrient storage that may be regulated by 825 unique environmental factors. Thus, we recognize multiple caveats within the trait-based model but expect that it improves our ability to link environmental and phytoplankton stoichiometry variation.

Conclusions
We find that processes that affect nutrient supply in oligotrophic gyres, such as the 830 strength of the thermohaline circulation, are particularly important in setting pCO2,atm but via a complex link with C:Pexport. By explicitly modeling the shallow overturning circulation, we showed that increased export in the tropics, which might be influenced by increased atmospheric iron dust deposition, may lead to increases, rather than decreases, in pCO2,atm. Increased [P] drawdown in the tropics shifts export away from the subtropical gyres, and 835 changes the mean export C:P in the low-latitude ocean. Additionally, we find that it is difficult to separate nutrient supply and temperature controls on marine phytoplankton stoichiometry, carbon export, and pCO2,atm and we need better physiological experiments and field data to fully understand the relative impact of the two factors. Nevertheless, it is likely that both play a key role in regulating phytoplankton stoichiometry, C:Pexport, and 840 ultimately ocean carbon cycling.
Author Contribution: ARM -creation and analysis of the box model and primary writer of manuscript. GIH -creation and development of the trait-based model, and writing. FWPassistance on the box model and editing of manuscript. SAL -assistance on the trait-based 845 model and editing of manuscript. ACM -assistance on both models and writing of manuscript.
Competing Interests: The authors declare that they have no conflict of interest.