Sensitivity of 21st-century projected ocean new production changes to idealized biogeochemical model structure

While there is agreement that global warming over the 21st century is likely to influence the biological pump, Earth system models (ESM) display significant divergence in their projections of future new production. This paper quantifies and interprets the sensitivity of projected changes in new production in an idealized global ocean-biogeochemistry model. The model includes two tracers that explicitly represent nutrient transport, lightand nutrient-limited nutrient uptake by the ecosystem (new production), and export via sinking organic particles. Globally, new production declines with warming due to 5 reduced surface nutrient availability, as expected. However, the magnitude, seasonality, and underlying dynamics of the nutrient uptake are sensitive to the light and nutrient dependencies of uptake, which we summarize in terms of a single biological timescale that is a linear combination of the partial derivatives of production with respect to light and nutrients. Although the relationships are non-linear, this biological timescale is correlated with several measures of biogeochemical function: shorter timescales are associated with greater global annual new production and higher nutrient utilization. Shorter timescales are also 10 associated with greater declines in global new production in a warmer climate and greater sensitivity to changes in nutrient than light. Future work is needed to characterize more complex ocean biogeochemical models in terms of similar timescale generalities to examine their climate change implications. Copyright statement. TEXT

prior work showing some parameter sensitivities are similar across models of varying complexity, results from our idealized approach may be widely applicable. In the course of our sensitivity study, we show that the magnitude of the new production response to climate change scales in proportion to a linear combination of the parameters quantifying the model's effective "biological timescale." Our study thus presents a new diagnostic, useful for studying physical-biological coupling in the context of a dynamic climate. This approach advances a conceptual framework via which inter-model differences in export-production 70 changes might be meaningfully deciphered. Our secondary motivation is to identify a model configuration suitable for studying process questions related to climate change at very high-resolution. Since computational costs scale in proportion to the number of simulated tracers, this amounts to finding a model capable of sufficient realism with the minimal number of tracers (Galbraith et al., 2015).
We describe our physical system, the development of the idealized model, and our analysis methods in section 3. In section 75 4, we describe the global rates, spatial patterns, and seasonal cycles of new production along with its controls and how these vary across parameter choices. This includes analyses of several regions that exemplify different physical climate perturbations.
Section 5 summarizes these findings and discusses the usefulness of this idealized model and the biological timescale.

Methods
We performed a set of idealized climate-change experiments with the Community Earth System Model (CESM) version 2.1 80 in an "ocean-sea-ice" configuration forced by atmospheric fields derived from reanalysis. In these experiments, the ocean and sea-ice models were integrated at a nominal 1 • × 1 • horizontal resolution; the ocean vertical grid included 60 layers, with 10 m resolution at the surface and 250 m at the ocean bottom of 5500 m. Surface forcing was applied as a prescribed atmospheric state, with a repeated annual cycle, based on the Coordinated Ocean-Ice Reference Experiment (CORE) Protocol (Large and Yeager, 2004).We added a pair of tracers representing idealized nutrient and phytoplankton, where production depends on 85 nutrient and light availability alone, not existing biomass, and plankton explicitly sink while being advected and mixed. The following subsections describe the details of the physical model runs (2.1), the development of the idealized biogeochemical tracers (2.2), followed by the method used for analyzing the causes of changes in production (2.3).

Timeslice experiments
To develop a process-oriented means of examining the response of new production to idealized changes in climate, we adopt 90 a "timeslice" approach. Rather than running a full transient integration, we perturb the model's initial conditions and the surface forcing to simulate a period representative of the future climate state. We thus run separate integrations designed to be representative of early-and late-century climate conditions. For the former, we begin from a state initialized from observations and integrate the model with forcing representative of a statistically normal annual cycle, i.e. a normal year (Large and Yeager, 2004); we perform a 20-year spin-up, which is sufficient to minimize interannual drift in the physical state, and then use 10 95 further years as our early-century timeslice. For our late-century timeslice, we adjust the initial ocean state and atmospheric forcing variables using anomalies computed from the fully-coupled CESM1 Large Ensemble (CESM-LE; Kay et al., 2015).
The CESM-LE includes 40 members integrated from 1920-2100; we use anomalies computed from the ensemble mean, quantifying the difference in ocean state and atmospheric forcing variables between 2000 and 2100. We then integrate with the normal-year forcing plus monthly ensemble-mean atmospheric anomalies from the CESM-LE, again using 10 years as our 100 timeslice. The LE has been examined in comparison to both CMIP5 (Alexander et al., 2018) and observations (Deser et al., 2017), with results showing similarity to future SST and observed climate variability, respectively. Using the century-scale mean changes from the CESM-LE allows us to represent the forced changes over the 21st century from the RCP8.5 scenario without the "noise" associated with natural inter-annual variability represented in any individual ensemble member.
Our resulting model runs have ocean temperature and salinity representative of early and late 21st century very similar to 105 the LE. Drift in these values within the decadal runs are small compared to either the imposed change between them or typical inter-annual variability in a coupled model. The atmospheric surface state has the same sub-seasonal variability in each year and both epochs and no inter-annual variability. Physical fields of interest for the biological impacts of climate change include changes in available light, vertical stratification, and ocean currents. We discuss only these fields that can impact our idealized biogeochemical model, leaving out other fields that do not have direct impacts, such as temperature or pH. The changes in these  1e). With this framework of the physical changes, we can consider their impacts on biological rates.

Model formulation
Our aim in developing a set of idealized tracers to represent new production and export is to have a minimal model which allows us to explicitly connect responses of these biological rates to the physical climate perturbation described above. This section describes the assumptions used to design the tracers and their mathematical form, followed by a sensitivity analysis in the 2000s climate and the choice of a limited set of parameters to analyze further. To explicitly represent the supply of  tracer can represent the phytoplankton that is created and follow it to depth. We assume an equivalence between new and export production due to the need for mass balance in steady state; thus, the annual supply of nutrient from depth, used in new production, is expected to match the downward flux of plankton and detritus (Eppley and Peterson, 1979).
In designing the nutrient tracer, we make three simplifying assumptions. First, we assume that the deep nutrient pool has 125 a fixed concentration, not dependent on explicit remineralization, which decouples the nutrient tracer from the export tracer.
Second, we assume that new production depends on the availability of this nutrient and light alone, not on the water temperature or on the existing plankton population that may be sustained by recycling of nutrients; this omits processes thought to be important in bloom-type events (Behrenfeld and Boss, 2014) but again keeps the nutrient and export tracers decoupled. Finally, we assume that the light available for new production in the mixed layer is an average of the light levels within the mixed 130 layer (as done in McGillicuddy Jr et al., 2003); below the mixed layer, productivity depends on the light at only the depth in question. This choice increases subsurface growth within the mixed layer, while allowing growth below the mixed layer depth.
With these assumptions, the reactions of the nutrient, N , are governed by the following equation: where µ 0 is the maximum growth rate (mmol N/m 3 day), Q is the nutrient limitation (nondimensional), L is the light limitation (nondimensional), k N is the half-saturation constant for the nutrient (mmol N/m 3 ), α is the sensitivity for the light limitation (m 2 /W), and S is the restoring source at depth. Light, I (W/m 2 ), decays exponentially from the surface value of the incoming short-wave radiation, but is averaged over the mixed layer depth (defined by maximum buoyancy frequency criterion, Large 140 et al. (1997)). Thus, light has a constant value in the mixed layer and an exponentially-decaying value below, which does allow growth below the mixed layer depth. For z ≤ −1000m, N is restored to a value of 20; this is our fixed nutrient pool. While the observed deep nitrate values vary on the basin-scale, this fixed N value at depth decreases drift in the solutions as only the upper ocean's nutrients must spin up. There is no flux through the air-sea or sea-land interfaces. The physical transport and mixing are done by the same mechanics as existing passive tracers in POP.

145
In designing the second tracer, we aim to explicitly represent the export from the surface to the deep ocean. This export is a combination of plankton and detritus, but should be the same total mass, on average, as the supply of nutrient upward. We do not differentiate between different types of sinking organic matter, and refer to them as a whole as particles. We assume that there is a constant sinking rate relative to the surrounding water and a small loss rate, which represents both remineralization into the recycled nutrient pool and removal through higher trophic levels. Particles, P , have their reactions governed according 150 to the following equation: where σ is the specific mortality rate of particles (1/days) and w s (m/day) is the vertical sinking rate of particles. There is no flux through the air-sea or sea-land interfaces. Again, advection and mixing are applied by the existing CESM-POP mechanics for passive tracers.

2000s Sensitivity
Our idealized tracers have five parameters: µ 0 , k N , α, σ, and w s . In order to identify reasonable values, we perform a sensitivity analysis of the first three under the early-century climate. The specific mortality rate and sinking rate of P are held fixed at w s = 0 and σ = 1/60days in this analysis, as they do not affect N or the production rate, and only modestly affect the horizontal spatial patterns of P : faster mortality reduces global-mean P , while faster sinking moves P deeper in the water To examine sensitivity to parameters and choose reasonable cases for the climate-change experiments, N is compared to near-surface World Ocean Atlas nitrate distributions (WOA N O 3 , Garcia et al. (2013)). P is compared to output from the CESM biogeochemistry model, BEC (Moore et al., 2013) run with the same model physics. Specifically, we compare nearsurface P to total plankton, the sum of its three phytoplankton classes in carbon units, converted to nitrogen units using a 170 stoichiometric ratio of 16/117. Although our particles represent both living matter and detritus, near the surface this P is most like newly-produced plankton which we expect to have the same spatial patterns as total phytoplankton. We choose this comparison to take advantage of the work done to fit BEC to observations for the same physics (see This relationship between parameters and near-surface N is the result of a balance between physical supply of N from depth and its consumption by production at a rate set by the parameters. If production is slower than supply, N will accumulate, 185 increasing the production rate and decreasing the vertical N gradient and therefore the supply rate until a near-equilibrium is reached. Thus, we can predict this relationship between parameter choice and near-surface N using the production function, The total production's dependence on N and I can be understood through the initial slopes of the productionnutrient and production-light curves, which are µ 0 /k N and µ 0 α: these describe how quickly production increases for an initial injection of each limitation; steeper initial slopes and faster initial production mean less accumulation of N before equilibration. 190 Thus nutrient will be higher for higher k N /µ 0 and 1/(µ 0 α); these functions describe a two-dimensional space where increases along either axis decrease productivity and increase near-surface nutrient concentration.
We define a biological timescale as a single term to summarize these two axes. First, we recognize that k N /µ 0 has units of days, so we normalize α, which has units W/m 2 , by dividing by to reach the same units as 1/k N , m 3 /mmolN, so that 1/(µ 0 α/α 0 ) also has units of days. We then add the results. The choice of this normalization factor is specific to this model, where the currency is nitrogen. Then we see that we have a biological timescale for new production based on the parameters chosen; we call this biological timescale τ bio : The format of this timescale reinforces the qualitative relationship between surface N and changes in µ 0 , α, and k N . We can 200 use τ bio as a metric of how quickly production might consume a new supply of nutrient; if we compared it to a physical timescale for nutrient supply rate, τ p , places where τ bio is shorter are likely to be nutrient-depleted and places where τ p is shorter are likely to be nutrient-replete. Across our 27 parameter cases, τ bio ranges from about 2 days to 2 years. In fig 2( τ bio is compared to global-mean surface N for the 27 cases and correlates well, r=0.96. This correlation is best for this and similar values of α 0 , e.g. 0.1 or 2 mmol N/WM, but is lower for e.g. 0.01 or 100 mmol N/Wm.

205
Meridional structures of P for these same integrations show more variation than N , which is clear in the more frequent intersections of the surface value curves and the variation in the latitudes of the maxima (figure 2c). The structure of the surface BEC phytoplankton has peaks near 60 • N, the equator, 45 • S, and along the Antarctic coast. The three cases with (α, µ 0 ) = (0.0125, 0.125) lack these peaks entirely, with a much smoother structure that we consider a poor match for observed behavior and will not be considered further. To determine the best of the remaining options, we measured the correlation 210 coefficient between the annual meridional structure of P and BEC phytoplankton, which ranged between 0.18 and 0.68 at the surface. The top half of these cases, 12 of the remaining 24 (see Table 1), were examined by eye for a good qualitative match at the surface, and will be used for further analysis. These cases have the correct general latitudinal structure of P and also have reasonable new (export) production rates between 5.5 and 7.6PgC/yr (figure 3a), which is within the 5-11PgC/yr rate in most literature (e.g. Cabre et al. 2015). 215 We will not be examining the explicit sinking of our particles. For those interested, consider sinking rates, w s , consistent with estimates for detritus: single cells sink at about 0.1m/day and aggregates can reach over 100m/day (Jackson and Burd, 1998). For decay rates, we suggest fixing these after choosing w s and adjusting them to have most production sink below the depth threshold of choice. For instance, with w s of 5m/day, σ =1/yr has 95% of annual production sink below 100m. We choose two cases for detailed analysis which have P structures capturing different aspects of the BEC surface phyto-220 plankton concentrations and quite different parameter and mean surface N values. The first case has surface P maxima near 45N and S and the equator, nicely matching those at the equator and 45S in BEC, but missing the 60N and Antartcic-coast peaks. This case has a small maximum growth rate (µ = 0.125mmol N/m 3 day), a moderate light sensitivity (α = 0.05), and a low nutrient threshold for growth (k N = 0.25). We consider this analogous to a small phytoplankton type or a phytoplankton community that has adapted to oligotrophic conditions, with plenty of light but low available nutrients. This case has a τ bio 225 of 162 days, corresponding to moderate mean surface N , slightly lower than WOA nitrate at the surface (see Appendix A for more details). We call this the 'slow' case.
The second case has surface P maxima near 60N, the equator, and the Antarctic coast, similar to those at 60N and the equator in BEC and capturing the increase toward the Antarctic coast, but missing the 45S maximum. In contrast to the first, this case has a fast maximum growth rate (µ = 2mmolN/m 3 day), a low light sensitivity (α = 0.2), and a moderate nutrient 230 threshold for growth (k N = 1). We consider this analogous to a large phytoplankton type or a phytoplankton community that has adapted to higher latitude conditions, with large seasonal cycles of light and nutrient availability. This case has a τ bio of 3 days, corresponding to a very low mean surface N , noticeably below observed nitrate concentrations. We choose this case for large contrast with the first and call it the 'fast' case. In the results section, we will focus on the changes in production under global warming for these two parameter cases of our idealized biogeochemical model and use the larger set of twelve cases for 235 context.

Light and Nutrient Controls on Production
Under global warming, changes in light and nutrient availability will vary both spatially and temporally. This section describes our method to untangle the effects of these two components on changes in new production. The form of production, R = µ 0 QL, allows for a decomposition and attribution of changes in productivity to changes in nutrient and/or light. As µ 0 , the maximum 240 growth rate, is held constant, we can examine simply the nutrient availability, Q, and the light availability, L, both of which are nondimensional and have values between zero and one. Using model output of the monthly-mean N and R, we compute the monthly-mean Q = N/(k N +N ) for each grid cell and the monthly-mean L = R/µ 0 Q. The change in QL can be decomposed as ∆QL = Q∆L + L∆Q + ∆Q∆L.
(6) 245 where ∆ indicates the change between the warmer, perturbed climate (2100) monthly values and those in the 2000s climate, and all terms can be averaged in space or time as needed. These difference terms are also nondimensional and have values between negative one and one. Due to the way we compute L from R and Q, there is no residual term-the above equation holds to numerical precision. We will use this decomposition to analyze both the spatial patterns of annual-mean production and seasonal cycles of production. When one of the three right-hand terms is most similar in size to ∆QL and highly correlated 250 in space or time, we will consider that term the main driver of changes in production.

Global statistics and dominant spatial patterns
For our twelve cases, global annual new production in the top 100m is 5.3-7.5 PgC in the early-century scenario and decreases to 4.8-6.2 PgC in the late-century scenario; production is higher for shorter τ bio in both epochs (fig 3ab). These patterns hold for our two exemplar cases, with the slow case having lower new production and a smaller decrease than the fast case. Our 12-case range of decreases in production, 9.5-19.5%, is similar to the 7-18% range of export decreases in CMIP5 (Fu et al., 2016) but larger decreases than most net primary production changes in CMIP5, -2 to -16% (Fu et al., 2016), or CMIP6, +6 to -12% .
Decreases in new production are partially determined by decreases in near-surface nutrient concentration. In all 12 cases, 260 global average N concentration in the top 100 m decreases by 15-22% in the late-century conditions (Fig. 3cd). Both initial concentrations and absolute reductions are smaller for shorter τ bio , but the reduction percentages highlight that the changes are somewhat insensitive to the varied light-and nutrient-limitation choices. These absolute reductions are 0.04-0.66, which is slightly smaller than the reductions of 0.66 ± 0.49 for CMIP5 RCP8.5 and 1.06 ± 0.45 for CMIP6 SSP5-8.5 . Despite the small absolute reductions in near-surface N concentration, the decrease in global new production 265 is larger for fast cases with shorter case timescale (Fig. 3b); these shorter τ bio cases have lower near-surface nutrient in the early-century conditions and higher slopes in their production-nutrient curve, making them more sensitive to changes.
The spatial patterns of annual production under early 21st century conditions demonstrate one effect of τ bio (fig 4a-b). The slow case has lower maximal values, as expected given its lower maximal growth rate, and smoother variations in space. The fast case has higher maximal values of new production, again set by its maximal growth rate, and sharper variations in space, 270 due to its greater need for nutrients in order to grow: production is more limited to locations with a fast physical nutrient supply.
Despite these differences, both cases match the general meridional patterns of observed phytoplankton (see e.g. Dasgupta et al. While we do not discuss spatial patterns in detail for the larger set of parameter cases, we do compute production by basin (supplement, Fig. B1). For every basin but the Arctic, production decreases for a warmer climate in all twelve cases, indicating that the qualitative results in most basins are not sensitive to the specifics of light and nutrient limitations. In the Arctic, which 285 has the least total production, the four cases with lowest τ bio have decreases in new production, while the remaining eight have increases. We will discuss this sensitive region further in Section 3.5.
As global reductions in production were matched by reduced nutrient concentrations, so too can spatial patterns in production be linked to nutrient concentrations and the vertical nutrient flux. Although the nutrient flux and new production do not necessarily align in space and time, we empirically find pattern correlations between the annual mean production and the 290 nutrient flux at 100 m as well as the average nutrient in the top 100 m (Fig. 5). The vertical flux is somewhat more noisy, hence the correlations of the absolute changes between present and future are only r=0.51, r=0.14 for the fast and slow cases, respectively. The average nutrient concentration in the top 100m reflects the integrated effects of the flux and productivity, hence the correlations are somewhat stronger. In the present-day, the pattern correlations between annual nutrient concentration in the top 100 m and production integrated to 100 m are r=0.56 and r=0.26 for the fast and slow biogeochemistry; for the future they 295 are r=0.58 and r=0.33, and for the percent changes r=0.88 and r=0.48, respectively. Across our 12 cases, the lower τ bio cases, which are more nutrient-limited, have stronger correlations between nutrient concentration in the top 100m and productivity.
These correlations are always stronger in the warmer climate, when all cases are more nutrient-limited.
While we have connected the production and its changes to nutrient concentrations, light effects are also important. Productivity is a product of functions representing the sensitivity of productivity to light L(I) and to nutrient Q(N) availability. As 300 described in the methods, this allows a decomposition of the changes in production into those caused by changes in nutrient availability, L∆Q, those caused by changes in light availability, Q∆L, and the covariance of the two, ∆Q∆L. Examining the spatial patterns of QL, its change, ∆(QL), and the components of that change (Fig. 6) provides more details on the controls of reduced production beyond the vertical nutrient supply and surface concentration. From the spatial fields, QL and ∆(QL) show the same spatial patterns as production and its changes (by definition, r = 1), although the latter is somewhat different 305 from the percent changes. From the components of the change, it is clear that L∆Q is the main control on production, r = 0.74 and r = 0.64 for fast and slow, which is consistent with our discussion of nutrient supply and concentration (Fig. 6bcgh). The change in light availability is a smaller contributor, r = 0.19 and r = 0.43 for fast and slow (Fig. 6di). Finally, the covariance is anticorrelated to the total change, r = −0.054 and r = −0.025 for fast and slow (Fig. 6ej). This component is largest near the equator, where it offsets increases in nutrient due to broader upwelling in the warmer climate.

310
From these analyses, we see that qualitatively, the spatial patterns of the response to the climate perturbation are similar across the biogeochemical cases under consideration, likely pointing to strong controls by the underlying physics via the nutrient supply. Quantitative differences exist; the spatial differences point to locations where results are sensitive to model parameters, while the global differences are consistent with the chosen parameters in that the faster cases, which have a higher nutrient utilization, have new production more correlated with the reductions of near-surface nutrient and its vertical supply.

Global Seasonal Cycle
Over much of the ocean, the seasonal cycle is the dominant mode of productivity variability due to seasonal variations in light and nutrients. We are interested in how the seasonal cycle will respond to our climate perturbation and whether this is sensitive to the case pararmeters.
We define the global seasonal cycle as monthly global spatial averages with a six-month offset between the northern and 320 southern hemispheres so that the Boreal and Austral seasonal cycles are in phase. Differences in the seasonal cycle of new production between the fast and slow cases are large and are a reflection of the different τ bio of the cases (Fig. 7b). Seasonal  flux. Productivity in the fast case quickly consumes the nutrients as they are supplied in the late winter and spring, whereas productivity in the slow case spreads out the nutrient consumption, achieving its maximum in mid-summer when the light is most plentiful (fig 7b).
Despite the qualitative differences between the present-day seasonal cycles, climate change has a qualitatively similar impact 330 on both the fast and slow cases (Fig 7a-c). The peak wintertime nutrient flux at 100 m (fig 7a) is reduced by 22% and 21% in the fast and slow cases, respectively. This reduction is consistent with, but may not be entirely caused by, reduced winter mixedlayer depths restricting the entrainment of nutrient-rich water from the thermocline (Fig 7f). Production is reduced during all months, but particularly in the latter half of the respective growing seasons. Thus, the growing seasons during the 2100s are modestly shorter with earlier peaks in both cases. It is true across all 12 cases that the largest reductions in new production are 335 after the peak new production. Thus, we can conclude that while the seasonal cycle can be very different across cases, for this model the global reduced production in a warmer climate follows a consistent pattern of a shortened growing season.
The causes of reduced production can be identified through the seasonal cycles of the components of ∆QL, all of which are small fractions of the maximum early-century QL values. In both slow and fast cases, winter and spring increases in light availability due to shallower MLD are offset by reductions in the nutrient availability and light-nutrient covariance, leading 340 to reduced production (Fig 7de). The minima of ∆QL are in the second half of the growing season in both cases, driving the shortened season. The slow case has negative values of all components at this time, with a largest reduction in nutrient availability (Fig 7d). By contrast, the fast case's minimum of ∆QL occurs while Q∆L is positive and is driven by reduced nutrient availability and especially the light-nutrient covariance (Fig 7e). This large effect of ∆Q∆L during the growing season is obscured in the annual mean ( fig. 6), where this factor is small over most of the ocean and spatially anticorrelated to the total 345 change, indicating the importance of considering the seasonal cycle for the mechanism of changes in production. In both our cases, reduced nutrient availability contributes to the shortened growing season, a consistent effect of the climate perturbation.

Regional Case Studies
We next examine three regions in detail: the downwelling South Pacific, the Arctic, and the Porcupine Abyssal Plain in the North Atlantic. These three regions exemplify the spatially-varied processes driving the climate response and the sensitivity of 350 that climate response across our model cases (figs 4,6). We emphasize that none of the model cases closely match observations in any particular region (see Appendix A for a comparison of N and observed nitrate), and that these examples are chosen for their qualitative differences. To contextualize these regions, we note their physical biomes; biome analysis has been used previously to identify ocean regions with similar biogeochemically-relevant physical characteristics across basins (Sarmiento et al., 2004;Cabré et al., 2015). A further discussion of biome analysis across the broader 12 parameter cases is available in 355 the appendix.
he downwelling South Pacific (red outline in fig 4), part of the permanently-stratified subtropical biome, is a region with small seasonal cycles, small climate-induced changes in physics, and broad losses in production across model cases. These losses of production are driven by the same processes in the fast and slow cases, suggesting an insensitivity to parameter choices. By contrast, the Arctic (north of 66.5 • N, the black line in fig 4), in the ice and marginal ice biomes, has large climate- induced physical changes, mainly in ice coverage. As noted earlier, the climate response in the Arctic is sensitive to parameter choice; this region is the most sensitive of those we discuss, with changes in annual new production having opposite signs in our two cases. The Porcupine Abyssal Plain (blue outline in fig 4), in the seasonally-stratified subtropical biome, has deep winter mixing in the present climate and large decreases in mixed-layer depths with the climate perturbation. This region is an example where the main driver (light or nutrients) of the climate-change induced reductions to productivity changes across 365 parameter space. These three regions provide illustrative examples demonstrating how physical climate changes create a range of possible responses in new production, as mediated by the parameters that modify the sensitivity of productivity to light and nutrient availability, and thus how sensitivity of projected changes in new production vary across the world oceans.

South Pacific
The downwelling South Pacific is defined by annual-mean downward vertical velocity at 100m depth in the South Pacific 370 between 10 and 35 • S. Maximum winter mixed layer depths reach only 120m, putting this region mainly in the permanentlystratified subtropical biome. Here, all biogeochemical parameter cases are qualitatively similar in both their current-climate seasonal behavior and their response to the climate perturbation, which we show is driven by reduced vertical nutrient supply.
Here, the seasonal cycles in both early-and late-century climates have high values of both upward nutrient flux and new production in the winter and spring (fig 8ab). Peak upward nutrient flux is in August in both cases in the early-and late-century 375 conditions. In the fast case, peak new production is in the same month, indicating a fast response to the supply of nutrients.
In the slow case, there is a one-month delay, indicating a slower response. These peak times do not change between earlyand late-century conditions; only the overall magnitude of rates decreases: in the slow case, annual production decreases by 25% from 0.11 to 0.086PgC/yr, while in the fast case, annual production decreases by 39% from 0.045 to 0.027PgC/yr. This is consistent with the global production being more reduced in the fast than slow case. This pattern holds across the 12 cases, 380 with all having peak production in August through October, no change in the peak timing with climate, and larger decreases in production for shorter τ bio .
The decrease in new production is due to the reduced availability of nutrients. The total change in production, ∆(QL), is nearly identical to the changes due to nutrient availability, L∆Q, with rms differences of 4 · 10 −5 and 2 · 10 −4 for fast and slow ( fig 8c); this near match holds across all 12 parameter cases. This is consistent with Cabré et al. (2015), who saw decreases 385 in production in low latitude nutrient-limited biomes across models. Our addition to that analysis is the direct quantitative connection to reduced nutrient availability that is not straightforward to compute for more complex biogeochemical models.
The reduced nutrient availability is due to reduced upward fluxes of nutrients, which we can exactly identify (fig 9ab). These fluxes are comprised of downward advective fluxes, upward fluxes from vertical mixing from KPP (Large et al., 1997), and upward fluxes from vertical effects of along-isopycnal mixing (Redi and GM). Both forms of parameterized mixing show 390 reduced annual fluxes in the warmer climate at 100m. Neither mixed layer depths nor effective vertical diffusivities change substantially, less than 10% for all MLD and less than 10% across 2/3 of the grid points in this region for diffusivity. Thus, reduced dN/dz is the driver of the decreased upward nutrient flux.
The changes in spatial-mean nutrient profiles (fig 9c) show a consistent pattern across the fast and slow cases: increased nutrient at 350-650m depth but decreases from 350m to the surface, including decreases near 100m. These mid-depth nutrient 395 changes are due to physical changes in circulation and/or mixing in the main thermocline. For diffusive transport, the decreases in N concentration above 350m mean that the vertical gradient of N near 100m, where we measure the fluxes, decreases. Thus, the decreased production in the downwelling South Pacific is being driven by decreased nutrients below the deepest winter mixed layers (about 120m), through decreased dN/dz causing decreased diffusive fluxes, decreased near-surface nutrient, and thereby decreased nutrient availability, Q. From these two examples it seems likely that the causes of reduced nutrient 400 availability in the warmer climate, which drives reduced production, are consistent across parameter cases.

Arctic
The Arctic region is defined by being within the Arctic circle (north of 66.5 • N), or having at least one day per year with no incoming solar radiation. This region has the largest sensitivity of changes in new production to model light and nutrient limitation, with the range across 12 parameter cases being a decrease of 17% up to an increase of 52%. Here, our two example 405 biogeochemistry cases have notably different responses to climate change, with the slow case having a 46% increase in annual production and the fast case having a 16% decrease.
Early-century seasonal cycles are similar for fast and slow cases. In this region, that is an upward flux of nutrients yearround and high new production rates in the summer, with a peak in May for the fast case and September for the slow case ( fig 10ab). Under late-century conditions, the peak production of the fast case is in the same month, with a slight increase in production earlier in the year and reductions later in the year which lead to a total decrease; in contrast, the peak production of the slow case is much higher and is earlier, now in July. These regional-mean production shifts are consistent with the different responses seen in the annual production maps (Fig 4), where the slow case has large increases while the fast case has moderate changes of both signs.
In the Arctic, increases in incoming light, with the seasonal maximum of the regional mean more than doubling (91 to 202 415 W/m 2 ), are due to reduced sea-ice (fig 1a) and have a large impact on changes in production. The mean increase in light availability in summer, quantified by Q∆L, its effect on changes in production, is nearly as large as the maximum QL in the 2000s for both cases (Fig 10cd). In the slow case, this increase in light allows for a large increase in production in the first half of the growing season, until nutrients are slightly depleted by that production, reducing production through lower nutrient availability in the second half (negative L∆Q and ∆Q∆L). This pattern is qualitatively consistent across slow cases with long 420 τ bio . In contrast, the increase in light availability in the fast case is offset by an associated reduction in nutrient availability, such that ∆Q∆L ≈ −Q∆L. Increased light in the spring leads to immediate increases in production (April), which uses enough nutrient to cause a dip in nutrient availability (May) before peak light availability (July), leading to decreased production in the summer and fall. This pattern is consistent across the four cases with shortest τ bio , which have decreased annual production in the warmer climate ( Fig B1). Thus, in the Arctic, the increases in light availability consistently increase production in spring, leading at some point to a biologically-driven decrease in near-surface nutrient. The resultant changes in production may take either sign and depend on the speed at which nutrients are removed, consistent with our τ bio .
The Arctic is a region where changes in new production under global warming is highly sensitive to the formulation of model production, both for our model and more complex ones. Vancoppenolle et al. (2013) found that CMIP5 projections of Arctic NPP were dependent on whether the Arctic reached a nutrient-limited state post sea-ice loss, with those models that reach a 430 nutrient-limited state having reduced production and others having increased production in the warmer climate. Our analysis provides a mechanistic hypothesis for these differences: shorter τ bio models have lower near-surface nutrients in early-century conditions, and therefore increases in production with increased light will more quickly lead to a nutrient-limited state.

North Atlantic
The Porcupine Abyssal Plain in the northeast North Atlantic is defined here to be 40-52N and 27-11W (cyan box in fig 4); 435 this includes the location of the Porcupine Abyssal Plain Sustained Observatory. This region is characterized by deep winter mixed layers and net downwelling. Under a warmer climate, the winter mixed layer depths are reduced; the spatial mean of the maximum decreases from 280 to 229m and the absolute maximum MLD is reduced from 669 to 590m (fig 11d). All 12 cases have decreased new production in the future, 9.55 to 26.75%, which is a smaller range than our other two regions but still larger than the global ocean; our two example cases have a decrease of 12% for slow and 21% for fast. In this region, as 440 for the Arctic, the mechanisms driving the changes vary.
Here, our two model cases show qualitatively different seasonal cycles with production mainly in the early spring for the fast case and spread over the summer for the slow case (fig 11b) similar to their respective global cycles. Both cases have little qualitative change between the early-and late-21st century. Most of the nutrient supply is in the winter for both cases (fig 11a), mainly due to wintertime entrainment associated with the mixed layer.

445
In the warmer climate, changes in production, ∆QL, depend on the relative impacts of increased light availability, Q∆L, in March and April and decreased nutrient availability, L∆Q, in all months; the signs of these drivers are consistent across cases, but the signs of the change in production are not. In the slow case, the changes in production with climate are an increase in spring and a decrease in summer and fall for a total decrease, closely following changes in light availability, 11c). The increased light in March and April corresponds to shallower mean and maximum mixed-layer 450 depths, respectively. Light availability decreases in the summer, due to increased mixed-layer depths in May and June. In the summer and fall, a higher portion of the reduction in production is due to reduced nutrient availability from both the lower winter peak in nutrient supply and the increased use during the spring production.
In the fast case, the same physical changes of the climate perturbation result in production decreases due to reduced nutrient availability in all months, ∆(QL) ≈ L∆Q (fig 11c). The winter increase in light availability has no noticeable impact, with 455 substantially lower winter-spring production co-occurring with the lower nutrient flux (fig 11a) and shallower monthly-mean mixed layer depth (fig 11d). This region's changes in new production are sensitive to parameter choices. Although total new production is reduced in both cases, reduced winter mixed layer depths can act to either increase or decrease spring production, depending on whether that production is more sensitive to light, as in the slower cases with longer τ bio , or nutrient availability, as in the faster cases with shorter τ bio . These differences highlight the usefulness of this idealized model, which allows us to 460 diagnose these drivers.

Conclusions
In order to study the sensitivity of the climate response of new production, we designed an idealized, two-tracer biogeochemical model that explicitly represents nutrient supply to the photic zone, new production, and the export of organic particles. The chosen simplifications for the production function allow for detailed analysis of the causes of changes in production but 465 eliminate or exaggerate other processes. First, dynamic phytoplankton concentration is omitted from the model; productivity is generally thought to scale with this concentration, so our seasonal results may be inaccurate in regions with strong blooms.
Second, the formulation of light limitation, using an average light concentration within the surface mixed layer, does not allow self-shading and enhances production near the bottom of the mixed layer. Finally, the simulations are not constrained to be realistic; the results provide information about the key parameter sensitivities.

470
From model integrations under early-and late-21 st century climate scenarios with 12 different parameter sets, we found that global production, near-surface nutrient concentrations, and projected changes in production were all connected to τ bio , derived from the initial slopes of the production-nutrient and production-light curves, which are the partial derivatives of production with respect to nutrient and light at the origin. A short τ bio indicates faster production and higher nutrient utilization, leading to  in a warmer climate. These reductions were linked to reduced near-surface nutrient availability and a shortening of the growing season in all cases. The percentage decreases in global new production were different by about a factor of two between the highest and lowest τ bio , a range similar to CMIP5. But, we find that the sign of the response (a reduction in productivity with warming) is the same for all of the light and nutrient limitation parameters that we considered.
We examined two exemplar cases, focusing on the similarities and differences in the seasonal, spatial, and regional responses 480 to climate change. Spatial patterns of changes in annual new production are similar between these two cases (pattern correlation r=0.26), with both being correlated to the changes in annual-mean nutrient availability (L∆Q), near-surface nutrient concentrations, and upward nutrient flux. Correlations were stronger for the fast case, which is more nutrient-limited. The ability to quantify the component of the reduction due to nutrient availability is unique to simple models like ours.
For more details on the drivers of the climate response and its sensitivity, we examined three regions. The South Pacific 485 region demonstrated the most consistent response to the climate perturbation. Here, changes in light and MLD are negligible.
Decreased nutrients in the upper thermocline drive lower vertical supply, lower nutrient availability, and lower production in all months for both exemplar cases. While there are still larger decreases for shorter τ bio , the mechanism is not sensitive to biogeochemical parameters.
From our detailed analysis of two higher-latitude regions, we found that compensation between changes in light and nutrient 490 availability have very different impacts for our runs with different production parameters. In the Porcupine Abyssal Plain, a reduced depth of the winter mixed layer acts to either increase or decrease spring production, depending on whether that production is more sensitive to light, as in the longer τ bio cases, or nutrient availability, as in the shorter τ bio cases. In the Arctic, larger increases in light availability due to sea-ice losses drive the largest sensitivity of new production's response to the climate perturbation, with different cases having opposite-signed annual-mean responses. Here, for fast cases (short τ bio ) 495 light-driven increased spring production reduces nutrient availability and thereby production later in the growing season to such an extent that annual totals are reduced. In contrast, slow cases can produce more throughout the growing season with little impact on nutrient availability. This analysis suggests a mechanism for the variation of projected Arctic production in CMIP5, where reduced production was associated with nutrient limitation (Vancoppenolle et al., 2013): if τ bio were diagnosed for those models, it may be that short τ bio cases have higher nutrient uptake, driving the nutrient limitation and reduced production.

500
Biological rates like τ bio are useful for explaining CMIP results more broadly. In those more-complex models, the effective τ bio varies in space because of changes in the phytoplankton community composition. It would be possible to quantify τ bio for these models empirically in a vertical column ocean through the individual injection of each nutrient at several concentrations to find the slope of the production-nutrient curve. Alternately, the new production rates from many locations in a global ocean configuration could be used to fit a production curve over all nutrients, with τ bio formed from the derivatives. Variations 505 between model results might then be related to their different τ bio , along with a comparison to the different physical rates.
We suggest that this reduced-complexity model and timeslice method may be suitable for high-resolution climate change process studies where computational cost is a limiting factor. Having only two tracers, this model is as inexpensive as possible while explicitly representing the supply of inorganic nutrient, new production, and sinking export. Given the range of produc-tion parameters that provide reasonable global results, one can choose the ones most suitable for the question of interest, such as approximately matching an ecosystem in a particular region or biome. While not considered here, the formulation of P allows for changes in export efficiencies under a warming climate to be studied through changing w s and σ.
Appendix A: Additional N-N O 3 comparison In our work, we were aiming to understand the sensitivity of the climate response and thus a wide range of possible model behaviors. To that end, our validation or comparison to observations and another, fitted, model, was quite simple. For further 515 context on how our two exemplar cases do (not) match observed behaviors, we show here the seasonal cycle of the 100maveraged N concentration compared to WOA nitrate concentrations. Figure  Plain region, where the slow case has N similar to nitrate but with a shift in the seasonal cycle. If future process work were to concentrate on an individual region, an analysis like this would allow for fitting of parameters.

Appendix B: Production in 12 cases
This appendix provides additional context of how our simple biogeochemical model's new production changes for different parameter values under the same physical climate perturbation. We provide a simple analysis of production changes across 525 ocean basins and biomes for the 12 parameter cases lightly considered in the main text. The fast case in the main text is here labeled by case timescale 3, the slow case by case timescale 162 (it is the rightmost 162 in figures). Global production rates and their changes in the warming climate were in figure 2.
First, we shown ocean basin averages of new production and its changes ( fig B1). These basins should be self-explanatory, but note that the Southern Ocean begins at 35 • S and the Arctic at 66.5 • N. For all 12 cases, the Southern Ocean is most 530 productive, likely related to our lack of iron limitation, and all basins except the Arctic show reduced new production in the warmer climate for all parameter cases. The largest percent losses are in the southern hemisphere basins, and parameter cases with shorter τ bio show larger reductions than slower.
Second, we show ocean biome averages of new production and its changes. Biome delineations are based on latitude, sea ice fraction, annual-mean vertical velocity at 100m, and maximum annual mixed layer depth. In the +/-5 • latitude band we have Upwelling regions in 5-30 • N 5-35 • S are the low latitude upwelling biome (LLU); above that, they are subpolar (SP) unless their ice fraction goes above 0.1 some month. Ice biomes are split for the northern and southern hemispheres, noted NI and SI respectively. Most production occurs in the ST SS, ST PS, and SP regions, with SP generally larger total than either ST region alone but smaller than their combination. All biomes' total production decreases in the warmer climate, by 3-60%. Production rates per area are larger in the Eq-U, ST-SS, SP, and SI regions than the Eq-D, ST-PS, LLU, and NI regions. Most biomes' production rate is lower in a warmer climate except Eq-D which increases for some parameter choices. From the changes in both production rate and annual production of each biome, it appears that the equatorial and ice regions' changes are most sensitive 545 to biogeochemical parameter choices, while the permanently stratified subtropics appear least sensitive. These sensitivities do not always follow the pattern from basin or regional analyses of faster timescale cases having larger reductions in production; an explanation is outside the scope of this analysis.
Biomes also shift in extent: Eq-U expands, decreasing the Eq-D area; ice, subpolar, and LLU contract (ice by over 20%,