A Bayesian model for assessing organic matter supply in complex marine food webs using amino acid stable isotope analysis

Model description

Food webs are complex networks in which organic matter is assimilated from basal resources and transferred across trophic levels in a series of predator–prey interactions. When considering the origin of the organic matter contained within an individual organism, various possible supply pathways with distinct sources can be envisioned. The OMSM described here uses a Bayesian framework designed to identify the basal source(s) of organic matter to an organism (or group of organisms) while also quantifying the number of trophic steps involved in the supply pathway. It achieves this by modeling the relationships between the isotopic compositions of amino acids measured in consumers (Y_j,k) and those measured in potential organic matter sources (X_j,i), where j represents a specific tracer (e.g., the δ¹⁵N value of an amino acid), i a possible organic matter source, and k an individual consumer sample. The model consists of four main components, each described in detail below (Fig. 1).

1. Designation of prior probability density functions (priors) for model parameters (e.g., mixing coefficients)

2. Incorporation of organic matter source data, X_j,i

3. Process model

   1. Linear mixing model

   2. Trophic discrimination model

4. Fitting of the model to consumer data, Y_j,k

Various parameters are specified or estimated within the model to represent the trophic structure of the food web, the relative contributions of organic matter sources, and the behavior of isotopic tracers within the system (Table 1). Markov Chain Monte Carlo simulation (MCMC) is used to explore the parameter space and identify the most probable parameter sets based on the observed data and defined prior information. Importantly, this approach allows the model to identify the breadth of possible and likely solutions, even when the mixing model is underdetermined and a single precise solution does not exist (Phillips & Gregg, 2003). The output of the MCMC is a set of posterior probability density functions, referred to hereafter as “posteriors,” which represent the most likely values for each parameter (Table 1).

Part 1 – Prior Information: Prior probability density functions are a necessary component of all Bayesian models and must be specified for all unknown parameters in the model. Priors are often defined to be intentionally vague, but they can also be used to bias the model towards or against specific solutions (e.g., specific mixing proportions) based on available knowledge that is independent of the main data inputs (e.g., the relative abundance of dietary sources in the environment or findings from stomach content analysis). While the default behavior of the OMSM is to use minimally informative priors (as described in Section ‘Choice of Tracers and Organic Matter Source Separation’), informative priors can be specified for mixing coefficients for each organic matter source (f_i) as in previous Bayesian mixing models (Moore & Semmens, 2008; Stock et al., 2018), as well as for the mean (μ_j,i) and standard deviation (σ_j,i) of each tracer within each organic matter source group. The construction of informative priors closely follows procedures used for MixSIAR and is thoroughly discussed by Stock et al. (2018). Unlike MixSIAR and other models, priors cannot be defined for trophic parameters (e.g., food web length, protistan trophic steps) due to the way that they are explicitly calculated within the model.

Part 2 – Organic Matter Source Data Incorporation: The model is provided with tracer values (here δ¹⁵N values of individuals amino acids) for samples representing each potential organic matter source. It assumes that the frequency of tracer values for each source can be described by a normal distribution characterized by a mean (μ_j,i) and standard deviation (σ_j,i). It is essential that sufficient data are available to accurately characterize the natural variability within each organic matter source. During the MCMC, the values for the mean (μ_j,i) and standard deviation (σ_j,i) of source tracer values are allowed to vary, approximating the variability contained in the organic matter source data. As a result, the model inherently incorporates natural variability in tracer values for each organic matter source without requiring explicit propagation of σ_j,i through subsequent equations. An example workflow for the selection of organic matter source data is presented in Section ‘Determining Possible Sources of Organic Matter’.

Part 3 – Process Model: The process model consists of two components: a linear isotope mixing model and a trophic discrimination model. Distinct sets of tracers are used to solve mixing equations (mixing tracers) and trophic equations (trophic tracers), requiring independent selection of the tracer sets best suited to quantify each process. To maintain the integrity of the model, there must be no overlap between mixing and trophic tracer sets.

The mixing model (3A in Fig. 1), which is similar to Bayesian mixing models implemented in other studies (Stock et al., 2018; Govan et al., 2023), describes how potential sources of organic matter combine to form the base of the food web. While this mixture of organic matter sources is often referred to as the “base” of the food web, we acknowledge that it may include detrital or non-algal material that does not strictly have a trophic position of 1. The model assumes linear mixing dynamics, such that the isotopic composition for any amino acid at the base of the food web can be estimated as: (1)${\displaystyle {\mathrm{\mu }}_{j,k,\text{base}}=\sum _{\mathrm{all}i}{f}_{i,k}{\mathrm{\mu }}_{j,i}=f_{1,k}{\mathrm{\mu }}_{j,1}+f_{2,k}{\mathrm{\mu }}_{j,2}+\dots }$

where the mixing coefficient, f_i,k, represents the fractional contribution of organic matter sources, i, to the base of the food web for an individual consumer, k, and μ_j,base is the average value of tracer, j, at the base of the food web. While mixing coefficients (here f_i,k) are the only unknowns in a typical mixing model, in the OMSM the tracer values at the base of the food web (i.e., “the mixture” or μ_j,base) are also unknown. The OMSM addresses this additional unknown by relating μ_j,base to tracer values measured in zooplankton via the trophic discrimination model.

The trophic discrimination model (3B) describes how the values of mixing tracers change, or are conserved, during trophic transfer. To implement this, trophic tracers are first used to estimate three key variables: the total food web length (FWL), the number of protistan trophic steps (PTS), and the number of metazoan trophic steps (MTS). Controlled feeding studies of zooplankton suggests that certain amino acids, such as the δ¹⁵N value of alanine (δ¹⁵N_Ala), exhibit a relatively constant rate of trophic discrimination across both metazoan and protozoan trophic steps (Décima et al., 2017). This discrimination is quantified using a trophic discrimination factor (Δ¹⁵N_AA), which represents the increase in δ¹⁵N_AA per trophic level (McClelland & Montoya, 2002; Chikaraishi et al., 2009; Bradley et al., 2015). By comparing the δ¹⁵N_Ala value measured in a consumer sample (δ¹⁵N_Ala,k) with the mean value at the base of the food web (μ_Ala,k,base), the total food web length for that sample can be estimated using the following equation: (2)${\displaystyle {\mathrm{\mu }}_{\mathrm{FWL},k}=\frac{{\mathrm{\delta }}^{15}{\mathrm{N}}_{\mathrm{Ala},k}-{\mathrm{\mu }}_{\mathrm{Ala},k,\text{base}}}{{\Delta }^{15}{\mathrm{N}}_{\mathrm{Ala}}}}$

where Δ¹⁵N_Ala is the trophic discrimination factor for δ¹⁵N_Ala. Alternatively, controlled feeding studies suggest that other amino acids (e.g., δ¹⁵N in glutamic acid/glutamate, or δ¹⁵N_Glx) only exhibit significant trophic discrimination during metazoan trophic steps (Gutiérrez-Rodríguez et al., 2014). Following the same logic used to estimate FWL, we can estimate metazoan trophic steps using the following equation: (3)${\displaystyle {\mathrm{\mu }}_{\mathrm{MTS},k}=\frac{{\mathrm{\delta }}^{15}{\mathrm{N}}_{\mathrm{Glx},k}-{\mathrm{\mu }}_{\mathrm{Glx},k,\text{base}}}{{\Delta }^{15}{\mathrm{N}}_{\mathrm{Glx},\text{meta}}}}$

where Δ¹⁵N_Glx,meta is the trophic discrimination factor for δ¹⁵N_Glx specific to metazoan trophic steps. Other tracers with comparable fractionation dynamics (i.e., those that exhibit constant trophic discrimination for estimating FWL or variable trophic discrimination across trophic groups for estimating MTS) can be substituted in these equations in place of δ¹⁵N_Ala or δ¹⁵N_Glx, provided the appropriate Δ¹⁵N_AA values are applied (see Table 2 and evaluation of TDFs below). An example workflow for the selection of Δ¹⁵N values are given in Section ‘Determination of δ¹⁵N_AA values’, while considerations and specific guidance for the selection of trophic tracers are given in Section ‘Choice of tracers and organic matter source separation’.

While uncertainty in μ_j,base is inherently incorporated into the model output, analytical uncertainty in δ¹⁵N_Ala/Glx,k and experimental uncertainty in trophic discrimination factors (Δ¹⁵N_Ala/Glx) must be explicitly propagated through Eqs. (2) and (3) to our estimates of food web length (FWL_k) and metazoan trophic steps (MTS_k). Using standard methods of error propagation, we can estimate the uncertainty in trophic parameters via the equations:

(4)${\displaystyle {\sigma }_{\mathrm{FWL},k}={\mathrm{\mu }}_{\mathrm{FWL},k}\cdot \sqrt{{\left(\frac{\sigma \left({\mathrm{\delta }}^{15}{\mathrm{N}}_{\mathrm{Ala},k}\right)}{{\mathrm{\delta }}^{15}{\mathrm{N}}_{\mathrm{Ala},k}}\right)}^2+{\left(\frac{\sigma \left({\Delta }^{15}{\mathrm{N}}_{\mathrm{Ala}}\right)}{{\Delta }^{15}{\mathrm{N}}_{\mathrm{Ala}}}\right)}^2}}$ (5)${\displaystyle {\sigma }_{\mathrm{MTS},k}={\mathrm{\mu }}_{\mathrm{MTS},k}\cdot \sqrt{{\left(\frac{\sigma \left({\mathrm{\delta }}^{15}{\mathrm{N}}_{\mathrm{Glx},k}\right)}{{\mathrm{\delta }}^{15}{\mathrm{N}}_{\mathrm{Glx},k}}\right)}^2+{\left(\frac{\sigma \left({\Delta }^{15}{\mathrm{N}}_{\mathrm{Glx}}\right)}{{\Delta }^{15}{\mathrm{N}}_{\mathrm{Glx}}}\right)}^2}}$

where σ represents the standard deviation associated with the indicated quantity. This uncertainty is incorporated into the model by drawing posterior values for FWL_k and MTS_k from normal distributions, with means defined by Eqs. (2) and (3), and standard deviations defined by Eqs. (4) and (5). The number of protistan trophic steps (PTS_k) is then calculated as the difference between FWL_k and MTS_k.

Now that we have used our trophic tracers to estimate PTS, MTS, and FWL, we can account for trophic discrimination in the mixing tracers between the base of the food web and consumers. For some tracers we can assume no trophic fractionation (termed conservative tracers; e.g., δ¹⁵N values of certain “source” amino acids (Popp et al., 2007) or, if used, δ¹³C values of essential amino acids (Fantle et al., 1999), yielding the simple equation: (6)${\displaystyle {\mathrm{\mu }}_{j,k}={\mathrm{\mu }}_{j,k,\text{base}}}$

where μ_j,k is the mean tracer value for some consumer sample. For δ¹⁵N_AA tracers where fractionation occurs identically during meta- and protozoan metabolism (termed constant trophic discrimination tracers; i.e., Δ¹⁵N_AA,meta = Δ¹⁵N_AA,proto), trophic discrimination can be accounted for uniformly throughout the food web via the equation: (7)${\displaystyle {\mathrm{\mu }}_{j,k}={\mathrm{\mu }}_{j,k,\text{base}}+{\mathrm{FWL}}_k\cdot {\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{meta}}.}$

Recognizing that, for some amino acids, isotope fractionation occurs differently during metazoan or protozoan metabolisms (termed variable trophic discrimination tracers; i.e., Δ¹⁵N_j,meta ≠ Δ¹⁵N_j,proto), fractionation occurring during meta- vs protozoan trophic steps is differentiated, yielding the equation: (8)${\displaystyle {\mathrm{\mu }}_{j,k}={\mathrm{\mu }}_{j,k,\text{base}}+{\mathrm{PTS}}_k\cdot {\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{proto}}+{\mathrm{MTS}}_k\cdot {\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{meta}}.}$

While variance in μ_j,base, PTS, MTS, and FWL affecting μ_j,k is inherently accounted for in μ_j,k due to the model’s structure, uncertainty in trophic discrimination factors (Δ¹⁵N_AA) still needs to be propagated through Eqs. (7) and (8). This uncertainty is then combined with analytical uncertainty in tracer measurements from consumer samples (hereon referred to as α_j,k). Using standard error propagation methods, we can estimate total uncertainty in the mean consumer tracer values (μ_j,k) with the following equations:

(9)${\displaystyle {\sigma }_{j,k}={\alpha }_{j,k}}$ (10)${\displaystyle {\sigma }_{j,k}=\sqrt{{\alpha }_{j,k}^2+{\left({\mathrm{FWL}}_k\cdot \sigma \left({\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{meta}}\right)\right)}^2}}$ (11)${\displaystyle {\sigma }_{j,k}=\sqrt{{\alpha }_{j,k}^2+{\left({\mathrm{PTS}}_k\cdot \sigma \left({\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{proto}}\right)\right)}^2+{\left({\mathrm{MTS}}_k\cdot \sigma \left({\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{meta}}\right)\right)}^2}}$

where Eqs. (9), (10) and (11) pertain to conservative (Eq. (6)), constant trophic discrimination (Eq. (7)), and variable trophic discrimination (Eq. (8)) tracers, respectively. Tracers are categorized into these three groups based on whether statistically significant differences between trophic discrimination factors in metazoans (Δ¹⁵N_AA,meta) and protozoans (Δ¹⁵N_AA,proto) have been observed in the literature (as in Section ‘Statistical Tests’), or if the tracer can be justifiably assumed to behave conservatively.

When the process model equations are combined by substituting Eq. (1) into Eq. (8), we obtain the full process model equation: (12)${\displaystyle {\mathrm{\mu }}_{i,k}=\sum _{\mathrm{all}i}\left[f_{i,k}{\mathrm{\mu }}_{j,i}\right]+{\mathrm{PTS}}_k\cdot {\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{proto}}+{\mathrm{MTS}}_k\cdot {\Delta }^{15}{\mathrm{N}}_{\mathrm{AA},\text{meta}}.}$

Inspecting Eq. (12) clarifies that the only true unknowns in the process model are the mixing coefficients, MTS, and PTS. This suggests that n sources + 1 independent tracers are required to fully constrain the model.

Part 4 – Model Fit: Results from model equations (i.e., the mean consumer tracer values, μ_j,k) are compared with the consumer data (Y_j,i) to assess how well the model parameters fit the observations. We assume that the tracer data in consumers (Y_j,i) follow a normal distribution with a mean (μ_i,k, Eq. (12)) and standard deviation (σ_AA,k, Eqs. (9)–(11)).

The output of the model is a set of posterior probability density functions (posteriors), which describe the most likely values for each model parameter (Table 1), given the data, the model equations, and priors. The model uses a Gibbs sampling protocol via the JAGS software package (Plummer, 2023) to identify the parameter values that produce the best model fit. Model posteriors are summarized by calculating the posterior mean and mode, as well as 50%, 75%, 90%, and 95% highest density intervals (HDIs), which represent the smallest interval that contains the respective percentage of posterior values. The mode, representing the peak of the posterior distribution, is the single most likely model solution, though all solutions within the 95% HDI range are considered reasonable and mathematically plausible solutions.

It is important to highlight the key differences between this model compared with other Bayesian mixing models. First, the number and nature of trophic levels comprising the food web (FWL, PTS, and MTS) is estimated within the model, whereas earlier models assume a single trophic transfer or can be forced to fit a defined number of trophic transfers. Second, this model explicitly differentiates between metazoan and protozoan trophic steps, allowing for each to be characterized by distinct Δ¹⁵N_AA values. This distinction is a critical advancement that enables more accurate modeling of trophic discrimination in the lower levels of planktonic food webs. While these features significantly expand the model’s applicability, they also increase the degrees of freedom, which can result in greater uncertainty and broader credible intervals in the model estimates.

OMSM implementation

The OMSM is implemented within a single R Markdown document that utilizes several functions defined in accompanying R scripts. The R Markdown guides the user through data importation, assessment of data properties, specification of model parameters, execution of the model, and analysis and visualization of the model output. The document and code are written in a general format, allowing them to be adapted to the specific application of other users. The OMSM itself is written in the BUGS (Thomas, 2006) and executed using the R implementation of the Bayesian analysis software JAGS (Plummer, 2023). All model code and dependencies are publicly available in a GitHub repository (https://github.com/CH-Shea/AA-CSIA-Organic-Matter-Supply-Model with DOI: https://doi.org/10.5281/zenodo.15548791).

The model requires three main data inputs. First, it uses tracer values measured in multiple samples from each potential organic matter source group. While the model can technically run with only two samples per group, three or more samples are preferred for more robust results. Second, it requires tracer values and associated measurement uncertainties for at least one consumer sample. If multiple consumer samples are included, the model fits separate mixing and trophic discrimination models for each sample. Finally, the model requires estimates for trophic discrimination factors for each tracer affected by trophic discrimination, along with their associated uncertainties. For tracers designated as having constant fractionation dynamics (following Eq. (7)), a single trophic discrimination factor is used (Δ¹⁵N_j,meta). For tracers with variable fractionation dynamics (following Eq. (8)), separate trophic discrimination factors are needed for protozoan (Δ¹⁵N_j,proto) and metazoan (Δ¹⁵N_j,meta) trophic steps. Recommended values for these parameters in marine planktonic ecosystems are given in Table 2.

OMSM interpretation

Because the outputs of the model are probability density functions (i.e., posteriors), the results must be interpreted probabilistically. Considering probabilistic intervals, such as highest density functions (e.g., HDIs), along with the distribution mode provides a convenient framework for hypothesis testing. The null hypothesis should be tailored to the specific research question. For example, if the question is: “Does surface POM contribute to mesopelagic organic matter supply?”, the corresponding null hypothesis would be “surface POM does not contribute to organic matter supply.” This hypothesis could be rejected with 95% confidence if the 95% HDI excludes zero. Regardless of this outcome, the lower and upper bounds of the HDI indicate the minimum and maximum proportional contribution of that source to the consumer food web, while the mode represents the most likely value. A similar logic applies to trophic parameters. For example, the HDI for PTS or MTS must exceed zero to confidently infer that these trophic steps are present in the food web. In this case, the HDI bounds provide the minimum and maximum estimated number of trophic steps, and the mode indicates the most likely value for each parameter.

Constraining sources of uncertainty to the OMSM

There are three main sources of uncertainty affecting the accuracy and precision of OMSM outputs:

1. Analytical uncertainty in consumer tracers values.

2. Observed uncertainty in organic matter source tracer values.

3. Observed uncertainty in tracer Δ¹⁵N values.

Although 2 and 3 have the greatest influence, the relative impact of each source depends on the structure of the food web being analyzed and can be assessed through careful consideration of Eq. (12).

Uncertainty in consumer tracers values (μ_i,k) will consistently affect model precision. However, analytical uncertainty in AA-CSIA measurements is generally small (standard deviation (SD) < 1‰). Moreover, because μ_i,k is not multiplied by any other term in Eq. (12), this uncertainty does not compound under different model conditions. As a result, this represents a relatively minor but consistent source of model uncertainty.

Variability in organic matter source tracer values can substantially affect model precision, particularly for the sources that contribute most to organic matter supply. While the multiplicative factor (f_i,k) associated with the organic matter source tracer value in Eq. (12) cannot exceed 1, within-source variability can be high (SD > 1‰), leading to significant error. Output precision can be improved by increasing sample size or reducing variance within source groups, especially for sources contributing most to organic matter supply.

Uncertainty in Δ¹⁵N values is another major contributor to model error, due both to high variability observed across controlled feeding experiments for many amino acids (Table 2) and to the fact that Δ¹⁵N values are multiplied by PTS, MTS, or FWL in Eq. (12). This means that Δ¹⁵N uncertainty has a greater effect as food web length increases. The most effective long-term mitigation will come from further research on amino acid fractionation in wild and laboratory organisms. In the meantime, users can reduce uncertainty by selecting tracers with lower observed variability in Δ¹⁵N values in the literature. Additional guidance on tracer selection and the effects of uncertainty in Δ¹⁵N values is provided in Section ‘Choice of Tracers and Organic Matter Source Separation’.

Methods & Data

Determining possible sources of organic matter

The first step in adapting the model is to identify all possible isotopically distinct sources of organic matter in the system. Previous studies at Station ALOHA (Gloeckler et al., 2018; Hannides et al., 2020), in the equatorial Pacific (Romero-Romero et al., 2020), and at Ocean Station Papa (Shea et al., 2023; Wojtal et al., 2023) have collectively identified three likely, isotopically distinct sources of organic matter that could support deep-sea food webs:

1. Surface POM – Particles collected by in-situ filtration above 50 to 100 m can be considered representative of fresh surface POM. This material may enter the mesopelagic food web directly by diel vertically migrating organisms that feed near the surface during the night, or indirectly through consumption of these migrators or their fecal pellets at depth.

2. Large mesopelagic particles – Particles >6 µm (Station Papa) or >51 µm (5 N, 8 N, Station ALOHA), as well as those collected in sediment traps below about 200 m depth, represent passively sinking particles in the mesopelagic zone. These particles are largely composed of fecal material or aggregated detritus.

3. Small mesopelagic particles – Particles 1–6 µm (Station Papa) or ∼1–51 µm (5 N, 8 N, Station ALOHA) collected below about 200 m represent the suspended or slow-sinking fraction in the mesopelagic zone. They are likely to have undergone some degree of microbial degradation and/or resynthesis.

Detailed descriptions of sample collection procedures and particle size classifications can be found in the primary literature (5N and 8N: Umhau et al. (2019); ALOHA: Hannides et al. (2020); Station Papa: Wojtal et al. (2023)).

AA-CSIA data from Station ALOHA clearly illustrates the general isotopic distinctions between these particle groups (Figs. 2 and 3). Small particles tend to exhibit the highest δ¹⁵N_AA values, while surface particles show the lowest, except in the case of δ¹⁵N_Thr values. Broader trends, criteria for sample selection and exclusion, and statistical comparisons among organic matter source groups are discussed in detail in the Supplemental Materials.

Determination of Δ¹⁵N_AA values

Trophic discrimination factors (Δ¹⁵N_AA) are used in the model to calculate trophic parameters (PTS, MTS, FWL) and account for trophic discrimination in mixing tracers. To determine an appropriate set of Δ¹⁵N_AA values for planktonic marine biota, Δ¹⁵N_AA values estimated in controlled feeding studies were aggregated (Gutiérrez-Rodríguez et al., 2014; McMahon & McCarthy, 2016; Décima et al., 2017) and filtered to only include aquatic organisms that excrete ammonia and have a TP of 2 or 3. In addition, one experiment on the fish Acanthopagrus butcheri fed vegetable meal was excluded for having anomalously high Δ¹⁵N_AA values. Experimental replicates were pooled such that the final compilation includes 21 experiments on four species of teleost, three species of marine gastropod, three species of pelagic crustacea (representing copepods and shrimp), the planktonic rotifer Brachionus plicatilis, and two species of planktonic protozoa (representing dinoflagellates and ciliates). Taxa were categorized as protozoa or metazoa, so that Δ¹⁵N_AA values could be compared between functional groups. The full data compilation is available in Table S2.

Several amino acids do not exhibit significant differences in Δ¹⁵N values between protozoan and metazoan taxa (Ala: p = 0.97, Pro: p = 0.21, Gly: p = 0.42, Ser: p = 0.21, Phe: p = 0.98, Lys: p = 0.29; Table 2). For these amino acids, Δ¹⁵N values are calculated as the mean across all experiments in the data compilation, including those of protozoans and metazoans. In the model, their δ¹⁵N_AA values are treated as tracers with uniform fractionation dynamics, following Eq. (7) for trophic discrimination.

In contrast, other amino acids show significant differences in Δ¹⁵N values between protozoan and metazoan taxa (Glx: p = 0.002, Asp: p = 0.02, Leu: p = 0.0001, Thr: p = 0.02; Table 2). For these amino acids, Δ¹⁵N values are calculated by averaging experiments conducted separately on protozoa and metazoa. In the model, their δ¹⁵N_AA values are treated as tracers with variable fractionation dynamics, following Eq. (8) for trophic discrimination. Regression analysis of δ¹⁵N_AA-Phe against δ¹⁵N_Glx-Phe values also yielded comparable Δ¹⁵N_AA values (Fig. S3), supporting their use in the model.

Two amino acids, isoleucine (Ile: p = 0.08) and valine (Val: p = 0.05), show marginal differences in their Δ¹⁵N values between protozoans and metazoans. Because of this uncertainty, we recommend excluding them from the model. Although they likely exhibit variable fractionation dynamics, their Δ¹⁵N_proto values are close to 0, resembling that of δ¹⁵N_Glx.

Δ¹⁵N_Thr exhibits the greatest variability across metazoan feeding experiments, ranging from −10.7‰to −1.4‰. In contrast, experiments with protozoans yielded more consistent values, ranging from −2.9‰to −1.5‰. This variability introduces considerable uncertainty when selecting an appropriate value for Δ¹⁵N_meta.

Given the large variability in Δ¹⁵N_AA values observed across controlled feeding experiments for many amino acids, we also used a linear regression approach similar to that of Bradley et al. (2015) to estimate Δ¹⁵N_AA values from wild zooplankton AA-CSIA data. This approach regresses δ¹⁵N_AA values against the known TP of the sample, such that Δ¹⁵N_AA values can be derived from the slope of that regression. While Bradley et al. (2015) based their TP estimates on stomach content analysis of teleosts, we instead used δ¹⁵N_Glx−Phe and δ¹⁵N_Ala−Phe as proxies for trophic position, assuming that our Δ¹⁵N values for these well-studied amino acids are reasonably accurate. Notably, this means we could not reliably estimate Δ¹⁵N_Ala or Δ¹⁵N_Glx using this method.

The regression analysis was based on a large compilation of AA-CSIA data measured in zooplankton and particles at equatorial Pacific sites at 5°N and 8°N (Romero-Romero et al., 2020), Station ALOHA (Hannides et al., 2020; Miller et al., 2025), and Station Papa (Shea et al., 2023; Wojtal et al., 2023). Particles collected above 100 m and zooplankton collected above 200 m were used to represent the surface-ocean food web, primarily composed of phytoplankton and their direct consumers. Phe-normalized δ¹⁵N_AA values (δ¹⁵N_AA-Phe) were used to minimize the effects of bacterially mediated isotope fractionation and account for isotopic baseline variability. Linear regressions were fit with δ¹⁵N_AA-Phe as the response variable and either δ¹⁵N_Glx-Phe or δ¹⁵N_Ala-Phe as the predictor. δ¹⁵N_Glx-Phe was used for amino acids that showed significant differences between their Δ¹⁵N values in metazoa and protozoa consumers in laboratory studies (Table 2), while δ¹⁵N_Ala-Phe was used for amino acids without such differences. Location was also included as a second, non-interacting predictor. Δ¹⁵N_AA was then calculated from the slope of the regression, m, using the equation: ${\displaystyle {\Delta }^{15}{\mathrm{N}}_{\mathrm{AA}}=\mathrm{m}\cdot \mathrm{TD}{\mathrm{F}}_{\mathrm{Glx}/\mathrm{Ala}}+{\Delta }^{15}{\mathrm{N}}_{\mathrm{Phe}}}$

where TDF_Glx/Ala is the Δ¹⁵N value of Glx or Ala relative to that of Phe. Additional details of these calculations and methods used for propagation of uncertainty are provided in the Supplemental Materials. The zooplankton and particle AA-CSIA data used, including all analytical uncertainties, and relevant code are available on GitHub (https://github.com/CH-Shea/Organic-Matter-Supply-Model/tree/main/Data) and published via Zenodo (DOI: 10.5281/zenodo.15548791).

Regression analysis comparing δ¹⁵N_AA-Phe against δ¹⁵N_Ala/Glx-Phe values in wild zooplankton produced Δ¹⁵N_AA values similar to those observed in controlled feeding studies (Fig. S3), providing confidence in the applicability of those values to wild zooplankton. Regression of δ¹⁵N_Thr−Phe against δ¹⁵N_Glx−Phe values produces a significant slope (ANOVA, F_2,89 = 103, p < 0.0001) with a reasonably strong fit (R² = 0.86). Because the resulting regression-based Δ¹⁵N_Thr value is in good agreement with the mean value observed in metazoan controlled feeding experiments but has less uncertainty (Table 2, Fig. S3), the model adopts a regression-based estimate of Δ¹⁵N_Thr = −5.9 ±1.5‰ for metazoans.

Choices for Δ¹⁵N_AA values used in the OMSM are summarized in Table 2.

Choice of tracers and organic matter source separation

Although the OMSM allows for the use of multiple amino acid tracers with non-zero TDFs, careful selection of the tracer set is essential and should be guided by the specific research questions. Notably, the considerations for selecting trophic tracers differ from those for mixing tracers.

For trophic tracers, two amino acids are needed: one with a constant trophic discrimination factor to solve Eq. (2) for food web length, and one with a variable trophic discrimination factor to solve Eq. (3) for metazoan trophic steps. Examination of Eqs. (4) and (5) suggests that choosing trophic tracers with the lowest relative uncertainty (std dev/mean) in their trophic discrimination factor minimizes the uncertainty of model outputs. Among amino acids with constant trophic discrimination factors (Ala, Pro, Gly, Ser, Phe, Lys), relative uncertainty is smallest for Pro (0.30) and Ala (0.40), making either one a strong candidate. We chose to use Pro because of its slightly lower relative uncertainty. For trophic amino acids with variable trophic discrimination factors (Glx, Asx, Leu, Thr), Glx has the smallest relative uncertainty (0.21), and was therefore selected as the most reliable variable TDF tracer. Accordingly, Pro and Glx were used as the trophic tracers in the model.

Selecting mixing tracers is more complex and requires evaluating multiple factors. First and foremost, a tracer’s ability to distinguish among different organic matter sources should be assessed, typically using Tukey’s pairwise comparison tests (or an alternative test if the data are not normally distributed). Effective mixing tracers should differentiate at least two, and ideally more, organic matter sources. Second, the degree of multicollinearity among potential tracers should be evaluated using correlation analysis (e.g., Pearson’s correlation coefficients, r) or by evaluating the effective dimensionality of the data set (Del Giudice, 2021). For n sources of organic matter, the best resolution is achieved using at least n-1 independent tracers. Therefore, identifying multiple, uncorrelated tracers is advantageous. However, if analytical error is a specific concern for a given tracer, including a second, redundant tracer may reduce the model sensitivity to analytical error, even if it does not provide additional information. Finally, inspection of Eqs. (10) and (11) reveals that uncertainty in the model output caused by fractionating tracers is influenced by the number of trophic steps (i.e., FWL, PTS, MTS) and by uncertainty in trophic discrimination factors. This implies that model uncertainty will be lower for shorter food webs, and that amino acid mixing tracers with lower absolute uncertainty in their Δ¹⁵N values should be prioritized.

For the Station ALOHA particle data, we find that the three organic matter sources (surface, large, and small particles) have distinct δ¹⁵N_AA values, which is visually apparent (Fig. 2) and statistically supported by permutational multivariate analysis of variance (PERMANOVA) (p < 0.001). Tukey’s pairwise comparison tests show that all 12 δ¹⁵N_AA values examined here differentiate at least two pairs of organic matter sources with >95% confidence, while 10 of them distinguish all three sources (see pairwise mean difference plots in Fig. S2). Despite this, visual inspection of the data shows that many amino acids display similar patterns in δ¹⁵N values across the dataset (Fig. 2). This is supported by correlation analysis (Fig. S3), which shows that many amino acid δ¹⁵N values are highly correlated, with a mean absolute Pearson’s correlation coefficient of |r| = 0.86. Trophic amino acid δ¹⁵N values are consistently highly correlated with one another (r = 0.87 - 0.99) and also correlated with source amino acid δ¹⁵N values (r = 0.85 –0.98). δ¹⁵N_Thr stands out as the least correlated to other amino acid δ¹⁵N values, with Pearson’s correlation coefficients ranging from 0.29 to 0.59, indicating it may provide uniquely independent information relative to other amino acids.

In addition to their ability to statistically distinguish all organic matter sources, δ¹⁵N_Phe and δ¹⁵N_Lys have the smallest associated uncertainties in their trophic discrimination factors (σ(Δ¹⁵N_Phe) = 0.5 and σ(Δ¹⁵N_Lys) = 1.2‰). Moreover, while highly correlated (r = 0.96), including both phenylalanine and lysine as mixing tracers could help buffer the model against analytical error arising from interference peaks that can occur around both amino acids during AA-CSIA (Fig. S4). These characteristics make them a logical pair of tracers to include in the model. While most other δ¹⁵N_AA values are also correlated with δ¹⁵N_Phe and δ¹⁵N_Lys values (Fig. S3), δ¹⁵N_Thr is notable for its weak correlation with these and other tracers, offering more orthogonal information. Principal component analysis of δ¹⁵N_Phe, δ¹⁵N_Lys and δ¹⁵N_Thr values (Fig. 4), reveals that 81% of the variance in the data is explained by the first principal component (PC1), which primarily reflects variation in δ¹⁵N_Phe and δ¹⁵N_Lys values and distinguishes small particles from surface and large particles. An additional 18% of the variance is explained by the second component (PC2), which primarily captures variation in δ¹⁵N_Thr values and separates surface and large particles. These three tracers have 1.7 effective dimensions, suggesting they can almost fully determine a 3-component mixture. Together, these findings support the combined use of δ¹⁵N_Phe, δ¹⁵N_Lys, and δ¹⁵N_Thr values to minimize model uncertainty and enhance the ability to resolve organic matter sources.

Finally, because phenylalanine has such a low trophic discrimination factor (Δ¹⁵N_Phe = 0.3 ±0.5‰), we propose that δ¹⁵N_Phe values can be treated as a conservative tracer in relatively short food webs, allowing us to avoid propagation of uncertainty in Δ¹⁵N_Phe, thus increasing the precision of model outputs. This assumption is supported by a complementary analysis presented in the Supplemental Materials, in which δ¹⁵N_Phe values were modeled non-conservatively. These analyses yielded similar, though slightly more variable, results (see Supplemental Materials at https://doi.org/10.5281/zenodo.15548791).

Procedures for model assessment using simulated zooplankton data

To assess model performance, the OMSM is tested using simulated zooplankton data, allowing model outputs to be compared with the known or “true” parameters used to generate the data. Particle data collected at Station ALOHA are used to represent realistic AA-CSIA data on organic matter sources to the mesopelagic zone. AA-CSIA data for 50 zooplankton samples (see Data S1) are then simulated by selecting realistic values for key ecological parameters (e.g., mixing coefficients, PTS, MTS, Δ¹⁵N_AA, etc.). Each of the 50 samples is assigned a unique set of mixing coefficients, representing a range of potential organic matter supply scenarios, drawn from a uniform Dirichlet distribution. These coefficients are applied to Eq. (1) to calculate the isotopic composition of each tracer at the base of the food web for each sample. Values for PTS and MTS for each sample are drawn from uniform distributions ranging from 0 to 1 and 1 to 2, respectively, yielding food web lengths ranging from 1 to 3. Using Eq. (12) and the mean Δ¹⁵N_AA values provided in Table 2, the δ¹⁵N_AA values for 50 simulated zooplankton samples are calculated. This data is then used as an input to the OMSM, enabling evaluation of model performance by comparing its output to the ecological parameters used to generate the simulated data. In diagnostic plots, the model output (or posterior distribution) for each model parameter (f_j,k, FWL, MTS, PTS, δ¹⁵N_AA,source, δ¹⁵N_AA,base, δ¹⁵N_AA,zoop) is compared to the “true” value selected or calculated during data simulation.

In model validation exercises we assume that we possess minimal prior knowledge. To do this, we use priors that will exert little influence over model outputs. For the mean organic matter source tracer value (μ_j,i) we use uniform prior distributions ranging from -100 to 100, for the standard deviation of organic matter source tracer values (σ_j,i) we use gamma prior distributions with shape and rate parameters set to 0.001 (Lemoine, 2019), and for mixing coefficients (f_i) we use a uniform Dirichlet prior distribution. When fitting the simulated zooplankton data, the OMSM is run using three parallel chains. Each chain includes a 5,000 step adaptation phase, a 10,000 step burn-in, and a 10,000 step sampling phase with a thinning factor of 10.

To assess the accuracy of the model’s estimated values for trophic and mixing parameters (posterior probability density functions), the known, simulated values for each parameter are compared to the mode of each posterior, identifying if the true value falls within the 50, 75, 90, or 95% HDI of the posterior. We also fit regressions of both the modelled value and the discrepancies against the true value to assess systematic error in model outputs.