Global patterns in marine organic matter stoichiometry driven by phytoplankton ecophysiology

We incorporated a macromolecular model of phytoplankton (CFM-Phyto) into the global ocean model (MITgcm). This combined model predicts cellular growth rate based on the macromolecular allocation, which in turn is used to determine the elemental stoichiometry of phytoplankton for the next model time step.

The phytoplankton component of the model is executed using the following algorithm, which is illustrated graphically in Extended Data Fig. 2: (1) relate the growth rate and elemental stoichiometry of phytoplankton based on the macromolecular allocation; (2) evaluate the possible growth rates under four different limiting nutrient assumptions and select the lowest rate: Liebig’s Law of the Minimum; (3) evaluate storage of non-limiting elements; (4) evaluate excess of non-limiting elements relative to maximum quotas; (5) based on that excess, evaluate effective nutrient uptake rate; and (6) evaluate the change in the elemental stoichiometry based on the balance between the growth rate and effective nutrient uptake rate. We describe the procedural details in the following text. Parameter values are listed in Extended Data Table 1. See ref. ²¹ for further details and justification of each equation in CFM-Phyto; here we repeat equations essential to explain the model used in the current study.

Connecting the elemental stoichiometry and the growth rate

The first step of the algorithm is to obtain the relationship between the current elemental stoichiometry and the growth rate (μ). To do that, we use CFM-Phyto²¹ (Extended Data Fig. 1). The model is based on the assumption of pseudo-steady state with respect to macromolecular allocation; in other words, the cellular-scale acclimation occurs rapidly relative to environmental changes. Laboratory studies show that macromolecular re-allocation occurs on the timescale of hours to days¹⁹. This is fast relative to the rates of environmental change in our coarse-resolution ocean simulations and so steady state solutions²¹ are used to relate growth rate, macromolecular allocation and elemental stoichiometry, as described in detail below. We first describe the case of N quota (here defined as Q_N; moles cellular N per mole cellular C) in detail, and then we briefly explain the case of P and C quotas as the overall procedures are similar. After that, we describe the case with Fe quota, which extends the previously published model²¹ for this study.

Relating N quota and growth rate

CFM-Phyto describes the allocation of N quota as follows, focusing on the quantitatively major molecules:

where Q_N is total N quota (per cellular C: mol N (mol C)⁻¹), the terms on the right-hand side are the contributions from protein, RNA, DNA, chlorophyll and N storage. We use empirically determined fixed elemental stoichiometry of macromolecules²¹ (Extended Data Table 1) to connect the macromolecular contributions of different elements (here C and P):

Here (Y_l^{j:k}) represents the imposed elemental ratio (elements j and k) for each macromolecular pool (l). (Q_{mathrm{C}}^x) and (Q_{mathrm{P}}^x) describe the contributions of macromolecule x to the total C quota (mol C (mol C)⁻¹) and P quota (mol P (mol C)⁻¹), respectively.

CFM-Phyto uses the following empirically supported relationship to describe (Q_{mathrm{P}}^{{mathrm{RNA}}}) (ref. ²¹):

where (A_{{mathrm{RNA}}}^{mathrm{P}}) is constant (below, A values represent constant except (A_{{mathrm{Chl}}}); see below), μ is growth rate (d⁻¹) and (Q_{{mathrm{P,min}}}^{{mathrm{RNA}}}) represents the minimum amount of RNA in phosphorus per cellular C (mol P (mol C)⁻¹). Substituting this equation into equation (3) gives:

In CFM-Phyto, we resolve three types of protein, photosynthetic, biosynthetic and other:

Photosynthetic proteins represent those in chloroplasts largely responsible for light harvesting and electron transport. The model assumes a constant composition of chloroplasts; thus, the amount of photosynthetic protein is proportional to the amount of chlorophyll²¹:

Biosynthetic proteins represent proteins related to producing new material such as proteins, carbohydrates, lipids, RNAs, DNAs and other molecules. The models use the following empirically derived relationship²¹:

Substituting equations (6)–(8) (in this order) into equation (5) leads to the following equation:

Empirically, chlorophyll depends on the growth rate and equation (10) accurately describes the relationship between the growth-rate dependences of chlorophyll under different light intensities²¹:

with (A_{{mathrm{Chl}}} = left( {1 + E} right)/v_I) and (B_{Chl} = m/v_I) with E (dimensionless) as a constant representing growth-rate-dependent respiration, and m (d⁻¹) describing maintenance respiration. v_I (mol C (mol C in Chl)⁻¹ d⁻¹) represents chlorophyll-specific photosynthesis rate based on an established function of light intensity I (μmol m⁻² s⁻¹)^21,57:

where (v_I^{{mathrm{max}}}) is the maximum chlorophyll-specific photosynthesis rate, e is the natural base and A_I is a combined coefficient for absorption cross-section and turnover time. Substitution of equation (10) into equation (9) leads to the following quadratic relationship between Q_N and μ:

where

Relating P quota and growth rate

Similarly, CFM-Phyto describes the relationship between the current P quota Q_P and μ. P is allocated to its major molecular reservoirs:

Similar to equation (7), with the assumption of the constant composition of photosynthetic apparatus, the model connects the amount of the chlorophyll to phosphorus in thylakoid membranes:

As for N allocation, substitution of equations (14), (4), (6), (7), (8) and (10) (in this order) into equation (13) leads to a quadratic relationship between Q_P and μ:

where

Relating C quota and growth rate

Similarly, CFM-Phyto describes C allocation as follows:

where Plip−Thy indicates P lipid in thylakoid membranes. The equation represents the allocation per total cellular C in mol C (mol C)⁻¹, so the sum of the macromolecules in C (Q_C) becomes 1. Using the imposed elemental ratios of macromolecular pools ((Y_l^{j:k})) we relate the elemental contributions:

Following the steps similar to those for the N and P allocations, substituting equations (14), (4), (6), (7), (8) and (10) (in this order) into equation (17) leads to the following quadratic relationship between cellular C quota Q_C (=1 mol C (mol C)⁻¹) and μ:

where

Relating Fe quota and growth rate

In order to capture global scale biogeochemical dynamics including the iron-limited high-nitrogen, low chlorophyll regimes, CFM-Phyto²¹ is extended to resolve Fe quota and allocation. The model is guided by a laboratory proteomic study⁵⁸ in which the major Fe allocations are to photosystems, storage and nitrogen-fixing enzymes (nitrogenase). As we do not resolve nitrogen-fixing organisms here, Fe allocation (mol Fe (mol C)⁻¹) represents only the first two:

As for equation (7) and equation (14), we relate the allocation of Fe to photosystems to the investment in chlorophyll, (Q_{mathrm{C}}^{{mathrm{Chl}}}):

This is a strong simplification because the pigment to photosystem ratio is observed to vary⁵⁹, but enables an explicit, mechanistically motivated representation of Fe limitation, which, a posteriori, results in global scale regimes of iron limitation that resemble those observed⁴³ (Extended Data Fig. 4). With equations (10), (19) and (20), we obtain the following relationship between Q_Fe and μ:

Evaluating the growth rate

We assume that the cellular growth rate is constrained by the most limiting element within the cell (and its associated functional macromolecules). Thus, at each time step and location, and for each cell type, the evaluation of growth rate is based on the following two steps: (1) computation of the growth rate for each element without storage; that is, the case when all of the elemental quotas are allocated to functional macromolecules; and (2) selection of the lowest growth rate among these; Liebig’s Law of the Minimum. For the first step, we define (mu _i) (i = C, N, P, Fe) as the growth rate, assuming that nutrient i is limiting. Under this condition, (Q_i^{{mathrm{Sto}}}) should be small as element i is allocated to other essential molecules. We assume that (Q_{mathrm{N}}^{{mathrm{Sto}}}) is also small under C limitation because N storage molecules are rich in carbon. With these assumptions, the solution for (mu _i) is obtained by solving the standard quadratic relationships of equations (12), (15) and (18) for N, P and C, respectively, neglecting any (Q_i^{{mathrm{Sto}}}) terms:

where Q_C = 1. For μ_Fe, equation (21) without (Q_{{mathrm{Fe}}}^{{mathrm{Sto}}}) leads to

Once the μ_i values are obtained, we determine the effective growth rate, μ, based on the lowest value, which identifies the limiting element based on current intracellular quotas:

Evaluating nutrient storage

In CFM-Phyto, non-limiting nutrients can be stored in an intracellular reserve²¹, reflecting commonly observed luxury uptake. Storage is evaluated as the difference between the total elemental quota (updated later) and the functionally allocated portion of that element:

Here (Q_i^{{mathrm{Non}}_{mathrm{Sto}}}) represents the contribution to element i by functional, non-storage molecules. For N, P and C, (Q_i^{{mathrm{Non}}_{mathrm{Sto}}}) is represented by the non-(Q_i^{{mathrm{Sto}}}) terms on the right-hand side in equations (12), (15) and (18), respectively:

Similarly, for Fe, from equation (21):

When there is N storage, (Q_{mathrm{C}}^{{mathrm{Sto}}}) must be recomputed to consider the allocation of C to it:

Evaluating the excess nutrient

Storage capacity for any element is finite and we define excess nutrient as a nutrient (N, P, Fe) that is in beyond an empirically informed, imposed maximum phytoplankton storage capacity. Excess nutrient is assumed to be excreted (see below). Excess of element i ((Q_i^{{mathrm{Exc}}})) is computed:

where (Q_i^{{mathrm{max}}}) is maximum capacity for nutrient i. For N, CFM-Phyto computes (Q_i^{{mathrm{max}}}) as a sum of non-storage molecules and prescribed maximum nutrient storing capacity according to model–data comparison²¹:

Laboratory studies suggest that when P is not limiting, the phosphorus quota maximizes to a value that is almost independent of growth rate^21,39,44. Storage of each element is finite and the upper limit to storage is imposed by specifying the maximum cellular quotas ((Q_{mathrm{P}}^{{mathrm{max}}}) (ref. ²¹) and also (Q_{{mathrm{Fe}}}^{{mathrm{max}}})) with size and taxonomic dependencies (for example, refs. ^27,41). Thus, the maximum storage is represented by the difference between the prescribed maximum quota and the actual quota²¹:

In the case where (Q_i^{{mathrm{Sto}}}) computed in the previous section exceeds (Q_i^{{mathrm{Sto}}_{mathrm{max}}}), the value of (Q_i^{{mathrm{Sto}}}) is replaced by (Q_i^{{mathrm{Sto}}_{mathrm{max}}}) and the difference is placed in the excess pool, (Q_i^{{mathrm{Exc}}}).

Computing effective nutrient uptake rate

One factor that influences the cellular elemental quota is the effective nutrient uptake rate (mol i (mol C)⁻¹ d⁻¹) of N, P and Fe, which we define as follows:

where V_i (mol i (mol C)⁻¹ d⁻¹) is nutrient uptake rate and the second term represents the exudation of the excess nutrient based on the timescale (tau _i^{{mathrm{Exu}}}) (d⁻¹). For V_i, we use Monod kinetics^60,61:

where (V_i^{{mathrm{max}}}) is maximum nutrient uptake, [i] (mmol m⁻³) is the environmental concentration of nutrient i and K_i (mmol m⁻³) is the half-saturation constant of i. Previous models have resolved the relationship between nutrient uptake and allocation to transporters^31,62. Here we do not explicitly resolve allocation to transporters, as proteomic studies indicate that it is a relatively minor component of the proteome compared with photosystems and biosynthesis in phytoplankton⁶³. Transporter proteins could be represented in a model with a finer-scale resolution of the proteome⁶⁴.

Differentiating small and large phytoplankton

In this model, ‘small’ phytoplankton broadly represent picocyanobacteria, which have high nutrient affinities and low maximum growth rates (for example, Prochlorococcus), whereas ‘large’ phytoplankton represent eukaryotes with higher maximum growth rates (for example, diatoms). The former are associated with a gleaner strategy adapted to oligotrophic regimes, while the latter are opportunistic, adapted to variable and nutrient-enriched regimes. To encapsulate this, the large phytoplankton have overall higher imposed (V_i^{{mathrm{max}}}) (~µ^maxQ_i), K_i and (v_I^{mathrm{max}}) than for the small phytoplankton (Extended Data Table 1), consistent with the previous models (for example, ref. ¹⁰). In addition, the larger cells are assigned a higher (Q_{mathrm{P}}^{{mathrm{max}}}) following the observed trends (Fig. 1 and Extended Data Table 1).

Computing the change in the elemental stoichiometry

The computation of the change in the elemental quotas is done based on the balance between the effective nutrient uptake rate and the loss of nutrient to the new cells:

This change in the elemental quotas based on the cellular processes and the passive transport of elements in phytoplankton by the flow field created by MITgcm governs the elemental stoichiometry of the next time step at a specific grid box, as in other versions of ecological models with MITgcm¹⁰.

Calculation of CV values

We computed the CV values based on the following equation:

where σ is the standard deviation and (bar x) is the mean. The purpose of this computation is to quantify the latitudinal variation of the averaged elemental stoichiometry. Thus, we used the averaged values for each latitude (as plotted in Fig. 2) for the calculation of σ and (bar x).

MITgcm-CFM

The biogeochemical and ecological component of the model resolves the cycling of C, P, N and Fe through inorganic, living, dissolved and particulate organic phases. The biogeochemical and biological tracers are transported and mixed by the MIT general circulation model (MITgcm)^35,36, constrained to be consistent with altimetric and hydrographic observations (the ECCO-GODAE state estimates)⁶⁵. This three-dimensional configuration has a coarse resolution (1° × 1° horizontally) and 23 depth levels ranging from 5 m at the surface to 5450 m at depth. The model was run for three years, and the results of the third year were analysed, by which time the modelled plankton distribution becomes quasi-stable. Equations for the biogeochemical processes are as described by equations and parameters in previous studies^10,38. Here, however, we include only nitrate for inorganic nitrogen, and do not resolve the silica cycle. We simulated eukaryotic and prokaryotic analogues of phytoplankton (as ‘large’ and ‘small’ phytoplankton). The eukaryotic analogue has a higher maximum C fixation rate for the same macromolecular composition and higher maximum nutrient uptake rates, but also has overall higher half-saturation constants for nutrient uptake. We used light absorption spectra of picoeukaryotes, which sits in-between small prokaryotes and large eukaryotes¹⁰. In MITgcm, the mortality of phytoplankton is represented by the sum of a linear term (m_l), representing sinking and maintenance losses, and quadratic terms representing the action of unresolved next-trophic levels^66,67, implicitly assuming that the higher-trophic-level biomass scales with that of its prey. We assumed that the latter term is small to avoid introducing additional uncertainties. Similarly, we do not resolve the stoichiometric effects of prey selection due to the nutritional status of prey, or viral partitioning of nutrients in the environment⁵⁰. Atmospheric iron deposition varies by orders of magnitude around the globe and has a large margin of uncertainty, including the bio-availability of the deposited iron, which in turn depends on the source and chemical history of the deposited material⁶⁸. To realize a realistic global net primary production, we doubled the atmospheric iron input from ref. ¹⁰; this factor is well within the uncertainty of the iron supply estimates. Each of the two phytoplankton groups has variable C:N:P:Fe as determined by the component macromolecules at each time step. The pools of C, N, P and Fe are tracked within the modelled three-dimensional flow fields.

Source: Ecology - nature.com