Bioconcentration as a key driver of Hg bioaccumulation in high-trophic-level fish

The ability of monomethylmercury (MMHg⁺) to bioaccumulate in seafood is of concern due to its neurotoxic properties. The challenge in understanding MMHg⁺ bioaccumulation lies in the fact that its levels in higher trophic levels result from both bioconcentration and biomagnification. Furthermore, Hg can occur in several chemical species, including Hg²⁺ and MMHg⁺, which both bioaccumulate. Although the dominant pathway for MMHg⁺ bioaccumulation into seafood is the bioconcentration of MMHg⁺ in primary producers and the subsequent biomagnification to higher trophic levels, other pathways can contribute to MMHg⁺ bioaccumulation. In this study, we use a fully coupled 1D water column Hg bioaccumulation model to quantify how total bioaccumulation of Hg²⁺ and uptake of MMHg⁺ from the water (bioconcentration) in consumers affects the bioaccumulation of MMHg⁺ in high-trophic-level fish. The study is performed in three setups representing hydrodynamic conditions representative of the North and Baltic Seas. We find that Hg²⁺ bioaccumulation does not influence the bioaccumulation of MMHg⁺ but the bioconcentration of MMHg⁺ plays an important role. Although direct bioconcentration accounts for < 15 % of MMHg⁺ bioaccumulation in cod, the cumulative effect of bioconcentration on all trophic levels increases the MMHg⁺ content of cod by 28 %–49 %. These results demonstrate that bioconcentration in consumers is essential to accurately model the bioaccumulation of MMHg⁺ at higher trophic levels.

1 Introduction

The element mercury (Hg) is presently included in the World Health Organization's list of the 10 substances of greatest concern (WHO, 2020). This is due to the capability of Hg to be methylated to form monomethyl mercury (MMHg⁺), a potent neurotoxin generated by microbial methylation of inorganic Hg. MMHg⁺ biomagnifies within aquatic food webs, accumulating in predatory fish to concentrations that can impair human neurological development upon consumption. For example, it is estimated that the consumption of MeHg contaminated seafood contributed to 61 800 premature deaths and caused economic damage of up to USD 2.87 trillion globally (Chen et al., 2025). This issue is expected to become even more significant as anthropogenic Hg emissions are projected to increase in the coming decades (Maria Brocza et al., 2024). Despite recent efforts, the bioaccumulation of MMHg⁺ in the marine environment is a complex topic and is not yet fully understood. Part of the complexity of understanding MMHg⁺ bioaccumulation and toxicity is that Hg can undergo speciation and occur in the environment in several chemical forms with distinct physical and chemical properties (Bieser et al., 2023). In particular, these include dissolved Hg (Hg²⁺), dissolved elemental gaseous Hg (Hg⁰), MMHg⁺, and dissolved dimethylmercury (DMHg).

Often both MMHg⁺ and DMHg are combined and are termed methylmercury (MeHg). The importance of DMHg is currently debated. Although DMHg is a common form of Hg in deeper water, there are no measurements in the North and Baltic Seas that would differentiate between DMHg and Hg⁰, and its role can therefore not be assessed in the model (Fitzgerald et al., 2007). However, given the rapid photodegradation of DMHg in natural water and that it is generally not assumed to bioaccumulate, DMHg is assumed not to significantly bioaccumulate in biota in the coastal area investigated in this study (West et al., 2022; Morel et al., 1998). The strong bioaccumulation of MMHg⁺, on the other hand, can be attributed to several reasons: MMHg⁺ can be absorbed by phytoplankton by cell-dependent factors, such as membrane channels (Garcia-Arevalo et al., 2024), and it can strongly bind to sulfhydryl (-SH) groups in organic material, notably cysteine, which traps toxic MMHg⁺ in the cell (Arnold et al., 1983). Additionally, it has been shown by Tesán-Onrubia et al. (2023) that plankton communities in the southern Mediterranean Sea have lower MMHg⁺ concentrations than plankton in the northern Mediterranean Sea; they linked this to changes in environmental conditions affecting bioconcentration.

Both Hg²⁺ and MMHg⁺ can bind to organic material and bioaccumulate in the marine food web (Mason et al., 1995). However, due to the higher toxicity and bioaccumulation potential of MMHg⁺, the bioaccumulation of MMHg⁺ is the most important concern and receives the most attention (Mason et al., 2012). The two main pathways in which MMHg⁺ content of biota can increase are bioconcentration and biomagnification.

1.1 Used terminology: bioaccumulation and bioconcentration, biomagnification

Bioaccumulation in the marine environment refers to the total increase in pollutants in biota compared to that in the water. This can be quantified in nature by measuring the concentration of pollutants in both water and biota and calculating the ratio. This is typically expressed as the bioaccumulation factor, BAF. For example, the bioaccumulation of MMHg⁺ in organisms i can be calculated based on observations as:

$$\begin{array}{}\text{(1)}& {\text{BAF}}_i^{{\mathrm{MMHg}}^{+}}={\displaystyle \frac{C_i^{{\mathrm{MMHg}}^{+}}}{C_{\mathrm{w}}^{{\mathrm{MMHg}}^{+}}}}\end{array}$$

In which,

$$\begin{array}{}\text{(2)}& {\displaystyle }\begin{array}{rl}& {\text{BAF}}_i^{{\text{MMHg}}^{+}}={\text{The bioaccumulation factor of MMHg}}^{+}\\ & \text{ for organism }i\left[\text{L}{\text{kg}}^{-\mathrm{1}}\right]\end{array}\text{(3)}& {\displaystyle }\begin{array}{rl}& C_i^{{\text{MMHg}}^{+}}={\text{The concentration of MMHg}}^{+}\\ & \text{ in organism }i\left[\text{ng Hg}{\text{kg}}^{-\mathrm{1}}\right]\end{array}\text{(4)}& {\displaystyle }\begin{array}{rl}& C_{\mathrm{w}}^{{\text{MMHg}}^{+}}={\text{The free concentration of MMHg}}^{+}\\ & \text{ in water}\left[\text{ng Hg}{\text{L}}^{-\mathrm{1}}\right]\end{array}\end{array}$$

Since the BAF can be based on field measurements, it is a commonly used metric to estimate the link between the concentrations of pollutants in seawater and those in biota. In this study, we are interested in separating the bioaccumulation into separate pathways: the direct uptake from the water (bioconcentration) and the increase in pollutants due to trophic interactions (biomagnification).

Bioconcentration, is the increase in the concentration of Hg in biota directly due to uptake from the water. Because the process of bioconcentration relies on the exchange of Hg between the dissolved phase and an organism, it depends on the surface area of the organic material that is in contact with the water. Because of this, small organisms, such as bacteria and phytoplankton, have a greater ability to bioconcentrate Hg (Mason et al., 1996; Pickhardt et al., 2006). However, the bioconcentration process is controlled by a variety of factors, and recent studies show that the bioconcentration of Hg²⁺ is constant when normalized for cell density, while the uptake of MMHg⁺ is affected by changes in cell density and biomass. This suggests that MMHg⁺ uptake is influenced by cell-dependent factors, such as the thickness of the phycosphere and the availability of transmembrane channels, while this is not the case for Hg²⁺ (Garcia-Arevalo et al., 2024). Bioconcentration is typically defined by the bioconcentration factor (BCF). The BCF for MMHg⁺ in organisms i can be calculated as

$$\begin{array}{}\text{(5)}& {\text{BCF}}_i^{{\mathrm{MMHg}}^{+}}={\displaystyle \frac{{\mathrm{BC}}_i^{{\mathrm{MMHg}}^{+}}}{C_{\mathrm{w}}^{{\mathrm{MMHg}}^{+}}}}\end{array}$$

In which,

$$\begin{array}{}\text{(6)}& {\displaystyle }\begin{array}{rl}& {\text{BCF}}_i^{{\text{MMHg}}^{+}}={\text{The bioconcentration factor of MMHg}}^{+}\\ & \text{ for organism }i\left[\text{L}{\text{kg}}^{-\mathrm{1}}\right]\end{array}\text{(7)}& {\displaystyle }\begin{array}{rl}& {\mathrm{BC}}_i^{{\text{MMHg}}^{+}}={\text{The concentration of MMHg}}^{+}\\ & \text{ in organism }i\text{ due to direct uptake from water}\\ & \left[\text{ng Hg}{\text{kg}}^{-\mathrm{1}}\right]\end{array}\text{(8)}& {\displaystyle }\begin{array}{rl}& C_{\mathrm{w}}^{{\text{MMHg}}^{+}}={\text{The free concentration of MMHg}}^{+}\\ & \text{ in the water}\left[\text{ng Hg}{\text{L}}^{-\mathrm{1}}\right]\end{array}\end{array}$$

Here, Hg could either refer to Hg²⁺ or MMHg⁺. Note that for consumers this would define the theoretical BCF. In nature it is typically only possible to measure the BCF in primary producers, as in consumers it would be impossible to separate between MMHg⁺ that is taken up directly from the water and MMHg⁺ that is ingested via food. Bioconcentration is the most important step in bioaccumulation and phytoplankton can have a BCF of MMHg⁺ between 2×10⁴ L kg⁻¹ and 6.4×10⁶ L kg⁻¹ (Gosnell and Mason, 2015).

Biomagnification is when MMHg⁺ reaches higher concentrations at progressively higher trophic levels. The biomagnification factor, the fractional increase in MMHg⁺ with each trophic level, is estimated to be 7.0 ± 4.9 (Harding et al., 2018; Lavoie et al., 2013). This means that in addition to the high concentration in MMHg⁺ in phytoplankton, there is a large increase in MMHg⁺ at every consecutive trophic level. Because of this, low trophic level animals such as copepods typically bioaccumulate less MMHg⁺ than higher trophic level animals such as predatory fish. Even though the direct uptake rate from the water by the copepod might be higher than that in fish, there are fewer trophic levels below the copepod that provide the opportunity to biomagnify MMHg⁺, typically resulting in lower MMHg⁺ concentration in the food of copepods than in the food of high trophic level fish, and thus a higher overall MMHg⁺ concentration in the higher trophic level fish than in the lower trophic level copepod. Many seafoods consist of high-trophic animals, such as cod, tuna, or marlin, which can have trophic levels between 4 and 4.8 (Nilsen et al., 2008; Sarà and Sarà, 2007). Biomagnification can increase the already high levels of MMHg⁺ in phytoplankton by up to another factor 11.9^4.8≈145 420. This is typically defined by the biomagnification factor, BMF, which can be calculated assuming steady state for organism i, preying on organism j for MMHg⁺ as:

$$\begin{array}{}\text{(9)}& {\text{BMF}}_{i,j}^{{\mathrm{MMHg}}^{+}}={\displaystyle \frac{C_i^{{\mathrm{MMHg}}^{+}}}{C_j^{{\mathrm{MMHg}}^{+}}}}\end{array}$$

In which,

$$\begin{array}{}\text{(10)}& {\displaystyle }\begin{array}{rl}& {\text{BMF}}_{ij}^{{\text{MMHg}}^{+}}=\text{The biomagnification factor for}\\ & \text{trophic consumption of organism }j\text{ by }i\left[\text{unitless}\right]\end{array}\text{(11)}& {\displaystyle }\begin{array}{rl}& C_j^{{\text{MMHg}}^{+}}={\text{The concentration of MMHg}}^{+}\\ & \text{ in organism }j\left[\text{ng Hg}{\text{kg}}^{-\mathrm{1}}\right]\end{array}\text{(12)}& {\displaystyle }\begin{array}{rl}& C_i^{{\text{MMHg}}^{+}}={\text{The concentration of MMHg}}^{+}\\ & \text{ in organism }i\left[\text{ng Hg}{\text{kg}}^{-\mathrm{1}}\right]\end{array}\end{array}$$

The biomagnification factor of MMHg⁺ is extremely high, Lavoie et al. (2013) estimates the diet‐weighted average BMF in marine samples for MMHg⁺ as 7.0 ± 4.9 while it is below 1 for Hg²⁺ in most cases (Seixas et al., 2014; Lavoie et al., 2013). This, combined with the higher toxicity of MMHg⁺ is the reason why the bioaccumulation of MMHg⁺ is of much higher concern than the bioaccumulation of Hg²⁺.

Overall the dominant pathway in the bioaccumulation of MMHg⁺ is the bioconcentration of MMHg⁺ in phytoplankton and consequent biomagnification. The importance of this route is quantified by Wu et al. (2019) using a meta-analysis. They find that the concentration of MMHg⁺ at the base of the food web predicts 63 % of the observed variability in high-trophic-level fish, while the remaining 37 % is controlled by factors such as the Dissolved Organic Matter (DOM) content and oligotrophy.

1.2 Current models

Multiple models have been developed to explain MMHg⁺ bioaccumulation in marine ecosystems. Key examples include trophic transfer (Schartup et al., 2018), base-level accumulation (Zhang et al., 2020), planktonic bioaccumulation in the Mediterranean Sea (Rosati et al., 2022), MMHg⁺ dynamics on the Beaufort Shelf (Li et al., 2022), and speciation and bioaccumulation in the North and Baltic Seas (Bieser et al., 2023).

In all of the previous models, bioconcentration of MMHg⁺ is included as it is an essential driver. It is concluded in Schartup et al. (2018) that the bioconcentration of MMHg⁺ in zooplankton contributes less than 15 % of total MMHg⁺ bioaccumulation. Consequently, in later models such as presented by Rosati et al. (2022) this interaction is not included because their model focuses on the base of the food web. The study performed by Li et al. (2022) includes the process of bioconcentration for invertebrates, but it is not included for vertebrates. This means that our model would be the first model to include bioconcentration at every trophic level.

The bioaccumulation of Hg²⁺ is much less studied and not incorporated in any of the above-mentioned models. This is because Hg²⁺ is much less toxic than MMHg⁺ and therefore comparably understudied. While data is limited, this raises the speculative question if the link between the bioaccumulation Hg²⁺ and MMHg⁺ is not underestimated as Hg²⁺ and MMHg⁺ are in active equilibrium in the water.

The ECOSMO-MERCY coupled system, which is used by Bieser et al. (2023) is the only coupled model that models the bioaccumulation of Hg²⁺ and MMHg⁺ at higher trophic levels such as fish, while incorporating bioconcentration at every trophic level. The version used by Amptmeijer et al. (2025) expands on this by adding a second fish of even higher trophic level. Because of this, the ECOSMO-MERCY coupled system, as described by Amptmeijer et al. (2025) is used in this analysis.

To clarify, the base case used in this study is identical to the base case of the 1D model presented in Amptmeijer et al. (2025), where the bioaccumulation of both Hg²⁺ and MMHg⁺ is represented through biomagnification and bioconcentration across all simulated trophic levels. The difference lies in the direction and scope of the analyses. In Amptmeijer et al. (2025), the model was evaluated for both carbon stocks and Hg²⁺ and MMHg⁺ dynamics to demonstrate that it reasonably reproduces the removal of Hg from the water column via bioaccumulation. This evaluation enabled an assessment of the feedback mechanisms of bioaccumulation on Hg cycling and the overall Hg budget. In contrast, the present study focuses on a sensitivity analysis of two specific model design choices introduced in Amptmeijer et al. (2025): the incorporation of consumer-level bioconcentration of MMHg⁺ and the bioaccumulation of Hg²⁺, and provides a targeted evaluation of their influence on MMHg⁺ bioaccumulation.

1.3 The hypotheses

While MMHg⁺ is more concerning than Hg²⁺ at higher trophic levels, Hg²⁺ can form up to 98 % of the bioaccumulated Hg in phytoplankton (Pickhardt and Fisher, 2007). This results in a large removal of Hg²⁺ during the phytoplankton bloom period. However, it is demonstrated by Amptmeijer et al. (2025), which analyzes the feedback of bioaccumulation on Hg cycling, that there is no change in average tHg and aqueous Hg caused by bioaccumulation, but that there is a seasonal variation in the aquatic tHg content due to bioaccumulation. This means that even if the average concentrations of tHg are not altered by bioaccumulation, there may still be an effect of Hg²⁺ bioaccumulation on MMHg⁺ bioaccumulation, as during the phytoplankton bloom tHg is reduced which could lead to a reduction of available MMHg⁺ for bioaccumulation. It could be theorized that as the ecosystem reduces tHg during the phytoplankton bloom, it would reduce dissolved MMHg⁺, as this is in active equilibrium with other Hg species and therefore reduce the availability of MMHg⁺ for bioaccumulation. Based on this, we propose our first hypothesis that the bioaccumulation of Hg²⁺ can lower the bioaccumulation of MMHg⁺ by removing Hg²⁺, which in turn cannot be methylated and accumulated as MMHg⁺.

The majority of MMHg⁺ present in higher trophic levels is derived from their dietary intake (Lavoie et al., 2013). It is often assumed that MMHg⁺ bioconcentration is not crucial for its bioaccumulation at higher trophic levels based on results such as those presented by Schartup et al. (2018). It is, for example, omitted from several Hg cycling and bioaccumulation models such as the model presented by Rosati et al. (2022), or not incorporated into higher trophic levels, as is the case in the model presented by Li et al. (2022). However, this assumption overlooks that bioconcentration occurs at all levels of the trophic chain. For example, if microzooplankton and mesozooplankton acquire 5 % of MMHg⁺ through bioconcentration, mesozooplankton will have 5 % less MMHg⁺ from its diet, which consists of microzooplankton, and another 5 % less due to the absence of bioconcentration, leading to a total reduction of 10 %. The second hypothesis is that MMHg⁺ bioconcentration in consumers significantly elevates MMHg⁺ levels at higher trophic levels. This concept has been previously suggested and studied by Wu et al. (2019). Their research found that the BCF in fish spans 3 to 7 orders of magnitude and greatly differs across studied sites; yet, they did find a strong correlation between BCF and MMHg⁺ concentration in fish. Thus, we are not the first to suggest that direct water uptake is a significant factor in MMHg⁺ bioaccumulation; rather, this study extends this understanding by quantifying the role of bioconcentration in all consumers on MMHg⁺ bioaccumulation in fish at higher trophic levels.

Studies that have analyzed the relative contribution of bioconcentration in the bioaccumulation of MMHg⁺ in fish found that in freshwater fine-scale dace (Phoxinus neogaeus), the bioconcentration accounts for up to 15 % of the total bioaccumulation of MMHg⁺ (Hall et al., 1997). A study by Wang and Wong (2003) found that in marine sweetlips (Plectorhinchus gibbosus) bioconcentration in fish can dominate total MMHg⁺ concentration if they eat food with low MMHg⁺ levels, while intake from food dominates total MMHg⁺ uptake when fish eat food with higher MMHg⁺ levels. This means that there will be an effect of MMHg⁺ bioaccumulation in higher trophic levels as it is a direct flux of MMHg⁺ into the organism. In this study, we want to expand this and quantify the cumulative effect of MMHg⁺ bioconcentration in all consumers. This allows us to evaluate if consumer bioconcentration is indeed a small percentage of total MMHg⁺, or if it is a major contributor to the total MMHg⁺ concentrations.

It is important to analyze these interactions using models as they cannot be fully understood using field and laboratory studies. This is because Hg²⁺ and MMHg⁺ are in active equilibrium and we cannot measure MMHg⁺ in a system where phytoplankton would not absorb Hg²⁺. The effect of bioconcentration of MMHg⁺ can also not easily be measured. It is possible to measure the direct uptake and release of MMHg⁺ by higher trophic levels from the water column. This is done, for example, to estimate the bioconcentration rates of Hg²⁺ and MMHg⁺ in sweetlips in the earlier mentioned study by Wang and Wong (2003). The challenge lies in that the source (bioconcentration or biomagnification) of MMHg⁺ cannot be measured in observational studies. Additionally, MMHg⁺ consumed by higher trophic levels is bioconcentrated in lower trophic level consumers, complicating the quantification of bioconcentration's full impact.

To test the two hypotheses that Hg²⁺ bioaccumulation decreases MMHg⁺ bioaccumulation and that MMHg⁺ bioconcentration in consumers significantly increases MMHg⁺ bioaccumulation in higher trophic levels, we quantify the effect of bioaccumulation of Hg²⁺ and the bioconcentration of MMHg⁺ on the bioaccumulation of MMHg⁺. We do this by running the fully coupled GOTM-ECOSMO-MERCY coupled system used in Amptmeijer et al. (2025) with and without the bioaccumulation of Hg²⁺ and the bioconcentration of MMHg⁺ in consumers. Then we analyze the bioaccumulation of MMHg⁺ at different trophic levels in these different scenarios and finally evaluate the importance of both interactions.

2 Methodology

To quantify the importance of Hg²⁺ bioaccumulation and MMHg⁺ bioconcentration on the bioaccumulation of MMHg⁺ we modeled the bioaccumulation of MMHg⁺ in three different scenarios using three idealized 1D water column models representing different hydrodynamic regimes typical for the North and Baltic Seas.

2.1 Modeled region

The first North Sea setup is the permanently mixed – Southern North Sea at ($\mathrm{54}\mathit{°}{\mathrm{15}}^{\prime }{\mathrm{00.0}}^{\prime \prime }$ N, $\mathrm{3}\mathit{°}{\mathrm{34}}^{\prime }{\mathrm{12.0}}^{\prime \prime }$ E). The 41.5 m deep location of this setup is characterized by having constant water-column mixing. This remixing of nutrients within the euphotic zone creates good conditions for phytoplankton growth. Additionally, since the water column is mixed during the bloom period, macrobenthos can feed directly from the phytoplankton bloom, which results in a high macrobenthos biomass. This results in a high biomass turnover rate and makes macrobenthos an important food source for predatory fish (Heip et al., 1992). The Southern North Sea is rich in nutrients, and the phytoplankton bloom can be light-limited. Diatoms typically dominate the phytoplankton bloom in spring until silicate limitations reduce their growth, and flagellates can become dominant (Emeis et al., 2015).

The second setup is the seasonally mixed – Northern North Sea at ($\mathrm{57}\mathit{°}{\mathrm{42}}^{\prime }{\mathrm{00.0}}^{\prime \prime }$ N, $\mathrm{2}\mathit{°}{\mathrm{42}}^{\prime }{\mathrm{00.0}}^{\prime \prime }$ E). This 110 m deep setup is only seasonally mixed. The Northern North Sea is still rich in nutrients, resulting in similar high phytoplankton growth, which is dominated by diatoms in spring and succeeded by flagellates in summer, similar to the Southern North Sea (Bresnan et al., 2009). A key difference between the Southern and Northern North Sea setups is that in the Northern North Sea setup, macrobenthos cannot feed directly on the bloom but predominantly feed on sinking detritus. This results in lower macrobenthos biomass and lower importance of macrobenthos in the diet of top predators. The reduced biomass of macrobenthos in the Northern North Sea is in line with observations (Heip et al., 1992).

The final setup is the permanently stratified – Gotland Deep ($\mathrm{57}\mathit{°}{\mathrm{18}}^{\prime }{\mathrm{00.0}}^{\prime \prime }$ N, $\mathrm{20}\mathit{°}{\mathrm{00}}^{\prime }{\mathrm{00.0}}^{\prime \prime }$ E). This setup is different from the two North Sea setups in several ways. First, the Baltic Sea, in general, is not limited by silicate, resulting in a dominance of diatoms in the phytoplankton bloom. In the Gotland Deep specifically, silicate limitation can occur, but diatoms will still be dominant throughout the bloom. Gotland Deep has a very low salinity (7 g L⁻¹), is strongly stratified, and can be eutrophied in phosphate. This results in good growth conditions for nitrogen-fixing cyanobacteria that can form a major part of the total phytoplankton biomass in autumn when nitrogen limitations limit the growth of other phytoplankton (Kahru and Elmgren, 2014). The presence of cyanobacteria can alter bioaccumulation because they can reduce dissolved Hg²⁺ to volatile Hg⁰, which increases Hg evaporation and therefore lowers the concentration (Kuss et al., 2015). This can reduce the average Hg content by up to 8 % (Amptmeijer et al., 2025). At the same time, the small size of the cyanobacteria gives them an extremely high surface to biomass ratio, resulting in a very high bioconcentration factor of MMHg⁺ (Pickhardt et al., 2006). Finally, the Gotland Deep has anoxic water below the thermocline; because of this, there is no macrobenthos (Conley et al., 2009).

Quantifying the importance of the bioaccumulation of Hg²⁺ and bioconcentration of MMHg⁺ in consumers on MMHg⁺ bioaccumulation into higher trophic levels under these idealized circumstances will provide a unique insight into the drivers of the bioaccumulation of MMHg⁺ and increase our fundamental understanding of this process. Additionally, it is important to quantify the importance of these interactions using lighter models because their inclusion in models comes at a cost. Especially, the implementation of the bioaccumulation of Hg²⁺ is done by adding one state variable to the model per biota functional group, or two state variables if the biomagnified and bioconcentrated Hg²⁺ is treated as separate variables, as is done in the model used in this study. While this is feasible without much concern in the 1D water column models that we use in this study, when running large earth system models, adding insignificant state variables becomes an unnecessary waste of computational resources which results in the wasteful expenditure of research funds and energy.

2.2 The GOTM-ECOSMO-MERCY model

Hypotheses are evaluated using the Generalized Ocean Turbulence Model (GOTM) (Bolding et al., 2021) that is coupled to the ECOSMO E2E ecosystem model (Daewel et al., 2019) and the MERCY v2.0 Hg speciation model (Bieser et al., 2023). The models are coupled using the Framework for Aquatic Biogeochemical modeling (FABM) (Bruggeman and Bolding, 2014). This setup is chosen because it has been used and evaluated in previous studies to analyze the bioaccumulation and cycling of Hg in the North and Baltic Seas, and it is the only fully coupled model to incorporate the bioaccumulation of Hg²⁺ and the bioconcentration of MMHg⁺ at higher trophic levels.

2.2.1 GOTM

GOTM is used to simulate the hydrodynamics of the 1D water column models. GOTM is a 1D turbulence model that computes the 1D version of the transport equation of temperature, momentum, and salinity. It does this while being nudged to observational datasets. GOTM simulations are designed using the iGOTM tool (https://igotm.bolding-bruggeman.com/, last access: 9 September 2019). This tool compiles the observational datasets used for the GOTM simulation and estimates the water depth based on the gridded bathymetry data ($\mathrm{1}/\mathrm{240}\mathit{°}$ resolution) (GEBCO Bathymetric Compilation Group, 2020), the ECMWF ERA5 data for meteorological data, the TPOX-9 atlas for tides ($\mathrm{1}/\mathrm{30}\mathit{°}$ resolution) (Egbert and Erofeeva, 2002), and for salinity and temperature, it uses the World Ocean Atlas (0.25° resolution) (Garcia et al., 2019). The state is solved every 60 s using forward Euler differential equations. The setups have 1 grid cell m⁻¹, and the variables are exported as daily means for the post-processing analyses.

2.2.2 ECOSMO E2E

The ecosystem model used in this study is the ECOSMO E2E ecosystem model (Daewel et al., 2019). The ECOSMO E2E ecosystem model is an intermediately complex ecosystem model that uses a functional group approach to estimate the biomass and carbon fluxes in the North and Baltic Seas. The version used here has three functional groups of phytoplankton: diatoms, flagellates, and cyanobacteria; two functional groups of zooplankton: microzooplankton and mesozooplankton; two functional groups of fish; and one group of macrobenthos. The basic version of the model is published by Daewel et al. (2019), but the version used here has some modifications to make it more suitable for bioaccumulation. This includes the addition of a second functional group of fish to represent high trophic level animals, the explicit resolution of the trophic level, and tuning of some of the carbon flux parameters such as growth rates and assimilation efficiencies. This is discussed in more detail in Amptmeijer et al. (2025).

2.2.3 MERCY v2.0

The MERCY v2.0 model links atmospheric Hg to MMHg⁺ in fish. It does this by estimating air-sea exchange and wet deposition of Hg based on the CMAQ-Hg model and calculating the marine cycling while taking into account marine speciation and bioaccumulation. It uses 35 state variables to estimate Hg speciation, transport, and bioaccumulation. The model estimates the partitioning of Hg²⁺ and MMHg⁺ into DOM and detritus, and bioaccumulation based on ecosystem parameters derived from the ECOSMO E2E ecosystem model (Bieser et al., 2023).

2.2.4 Bioaccumulation in the model

In addition to the bioconcentration, all consumer functional groups can take up Hg from the consumption of contaminated food. The uptake of Hg²⁺ or MMHg⁺ from food depends on the assimilation efficiency of the food and Hg species. After Hg has been assimilated from food, MMHg⁺ is released based on the mortality and respiration rate within the functional group, while there is an additional release rate for Hg²⁺. Since fish have a temperature-dependent respiration rate in the ECOSMO E2E model, this means that fish lose Hg from both bioconcentration and biomagnification faster in warmer water as their respiration, and thus Hg release rate, increases with temperature. The bioconcentration rates for zooplankton are based on Tsui and Wang (2004), and those for fish are based on Wang and Wong (2003).

The implementation of bioaccumulation is discussed and validated in more detail in Amptmeijer et al. (2025), but the core equations are discussed here as well for clarity.

The increase in bioconcentrated pollutant (Hg²⁺ or MMHg⁺) per day for a functional group is calculated based on the biomass concentration of the group, the uptake rate, and the concentration of the pollutant, while it is reduced with a rate that is the sum of the release rate of the pollutant and the loss of biomass from group g, from both biological loss (respiration and mortality) and predation. The change in pollutant p due to bioaccumulation can then be calculated using the following equation:

$$\begin{array}{}\text{(13)}& {\displaystyle \frac{\mathrm{d}{C}_{g,p}^{\mathrm{BC}}}{\mathrm{d}t}}=b_g\cdot C_p^{\mathrm{env}}\cdot r_{g,p}^{\mathrm{bc}}-C_{g,p}^{\mathrm{BC}}\cdot \left(r_{g,p}^{\mathrm{rel}}+r_g^{\mathrm{bl}}+\sum _{z=\mathrm{1}}^{n_z}{r}_{z,g}^{\mathrm{pred}}\right)\end{array}$$

• $C_{g,p}^{\mathrm{BC}}$ = Bioconcentrated pollutant p in group g [ng Hg m⁻³].

• b_g = Biomass of functional group g [mgC m⁻³].

• $C_p^{\mathrm{env}}$ = Environmental concentration of pollutant p [ng Hg m⁻³].

• $r_{g,p}^{\mathrm{bc}}$ = Bioconcentration rate for group g and pollutant p [m ³ mgC⁻¹ d⁻¹].

• $r_{g,p}^{\mathrm{rel}}$ = Release rate of pollutant p from group g [d⁻¹].

• $r_g^{\mathrm{bl}}$ = Biological loss rate for group g (mortality, respiration) [d⁻¹].

• $r_{z,g}^{\mathrm{pred}}$ = Loss of group g due to predation by consumer z.

• n_z = Number of consumer groups feeding on group g.

• z = Index for consumer groups (predators) of g.

• t = Time [d].

While the change in pollutant p due to biomagnification is also dependent on the predation and concentration of pollutants from both bioconcentration and biomagnification in the prey. Additionally, pollutant p is released via the turnover rate rather than the release rate, as is the case for bioconcentration. The change in pollutant p due to biomagnification can then be calculated as follows:

$$\begin{array}{}\text{(14)}& \begin{array}{rl}{\displaystyle \frac{\mathrm{d}{C}_{g,p}^{\mathrm{BM}}}{\mathrm{d}t}}& =\sum _{s=\mathrm{1}}^{n_s}\left(r_{g,s}^{\mathrm{pred}}\cdot a_{s,p}\cdot (C_{s,p}^{\mathrm{BC}}+C_{s,p}^{\mathrm{BM}})\right)\\ & -C_{g,p}^{\mathrm{BM}}\cdot \left(r_{p,g}^{\mathrm{to}}+r_g^{\mathrm{bl}}+\sum _{z=\mathrm{1}}^{n_z}{r}_{z,g}^{\mathrm{pred}}\right)\end{array}\end{array}$$

• $C_{g,p}^{\mathrm{BM}}$ = Pollutant p concentration in group g from biomagnification [ng Hg m⁻³].

• n_s = Number of prey groups consumed by g.

• s = Index for prey functional groups of g.

• $r_{g,s}^{\mathrm{pred}}$ = Predation rate of group g on prey group s [d⁻¹].

• a_s,p = Assimilation efficiency of pollutant p from prey s [unitless].

• $C_{s,p}^{\mathrm{BC}}$ = Pollutant p concentration in group s from bioconcentration [ng Hg m⁻³].

• $C_{s,p}^{\mathrm{BM}}$ = Pollutant p concentration in group s from biomagnification [ng Hg m⁻³].

• $r_{p,g}^{\mathrm{to}}$ = Turnover rate of pollutant p in group g [d⁻¹].

• $r_g^{\mathrm{bl}}$ = Biological loss rate for group g [d⁻¹].

• $r_{z,g}^{\mathrm{pred}}$ = Predation rate of predator z on group g [d⁻¹].

So the total concentration of pollutant p in ng Hg m⁻³ is:

$$\begin{array}{}\text{(15)}& C_{(g,p)}=C_{(g,p)}^{\mathrm{BC}}+C_{(g,p)}^{\mathrm{BM}}\end{array}$$

Since this tracks the pollutants per volume of water, the total bioaccumulation per biomass in ng Hg mgC⁻¹ is then calculated as

$$\begin{array}{}\text{(16)}& C_{(g,p)}^{\mathrm{bg}}={\displaystyle \frac{C_{(g,p)}}{b_g}}\end{array}$$

2.3 Performance of the GOTM-ECOSMO-MERCY model

The model performance is discussed in more detail in Amptmeijer et al. (2025), but the key metrics are summarized below. The model is generally consistent with observational data and the previously validated 3D ECOSMO E2E model in terms of biomass. Minor exceptions are that the chlorophyll a concentration in the Gotland Deep matches the Northern instead of the Central Baltic Sea, and that the fish biomass in the Gotland Deep is overestimated by 7 % compared to Thurow (1997). The model also predicts tHg content in phytoplankton, zooplankton, and fish 1 accurately, and the MMHg⁺ bioaccumulation in fish corresponds well with trophic interactions. A deviation is seen in the trophic level fish 2, which has a trophic level between 3.5 and 3.7 in the model, below the expected level for Atlantic Cod (between 4.0 and 4.2). Nonetheless, this level remains high, making fish 2 representative of high-trophic-level animals. The MMHg⁺ bioaccumulation in fish 2 is consistent with the observed bioaccumulation for its trophic level. Thus, with the above-discussed minor exceptions, the model simulates biomass, Hg speciation, and bioaccumulation in line with observations.

2.4 Post-processing analysis

The post-processing analysis is performed in R v4.4.1. Plots are generated using ggplot v3.5.0., and linear regression and statistics are calculated using ggpubr v0.6.0. The hypotheses are tested using two different tests. First, a Wilcoxon signed-rank test is performed using the wilcox.test function from the stats (v4.4.2) package in base R. This is a non-parametric statistical test that determines if there is a significant difference between the base case and the scenarios. We accept a p-value of < 0.05 to indicate that the scenario has a significant deviation from the base case. The scenarios are compared using a Bayesian t-test implemented via the ttestBF function from the BayesFactor package (v0.9.12-4.7). The Bayesian t-test computes the Bayes factor BF₁₀, defined as the ratio of the marginal likelihood of the data under the alternative hypothesis (different means) to the marginal likelihood under the null hypothesis (equal means). The analysis is performed using Jeffreys-Zellner-Siow priors (Zellner and Siow, 1980). A BF₁₀<1 indicates evidence for equal means, whereas BF₁₀>1 indicates evidence for a difference in the mean. Typically, a BF₁₀<0.1 or BF₁₀>10 is interpreted as strong evidence, and BF₁₀<0.01 or BF₁₀>100 as very strong evidence for equal or different means, respectively (Doll and Jacquemin, 2019).

2.5 Scenarios

The model is run in three different scenarios. The “Base case”, “No Hg²⁺ bioaccumulation” and “No MMHg⁺ bioconcentration”. The base case scenario is the same as the base case used in Amptmeijer et al. (2025). For the “No Hg bioaccumulation” setup, all uptake rates of Hg²⁺ are set to zero. For the “No MMHg⁺ bioconcentration” scenario, all consumer bioconcentrations of MMHg⁺ and all Hg²⁺ uptake rates are set to zero.

2.6 Sensitivity analyses

In order to further investigate how bioconcentration in consumers affects bioaccumulation of MMHg⁺, we performed a sensitivity analysis on the key drivers: the bioconcentration rate of consumers and the bioaccumulation rate of producers. To this extent, two sensitivity studies are performed. In the first sensitivity study, the bioconcentration rate in all consumers is multiplied by a scaling factor that is between 0.2 and 2.0 with 0.2 intervals. The effect of this on the bioaccumulation in fish 2 for the Gotland Deep is shown to visualize the impact. Then the relative contribution of bioconcentration in consumers on the bioaccumulation of MMHg⁺ in fish 2 is shown for all three setups. For the second sensitivity study, the same approach is used but the bioconcentration rate of producers is multiplied by a scaling factor. The effect of the consumer and producer bioconcentration scaling factor on the difference caused by consumer-level bioconcentration on total bioaccumulation is visualised. The data for both instances are plotted with an optimal fit for each of the three stations. The consumer bioconcentration factor is shown with a best fit assuming an asymptotic growth model $y=a\times (\mathrm{1}-e^{-bx})$, which is fitted using the nls function from the R stats package (version 4.4.2). For producer bioconcentration, a best fit with an asymptotic exponential decay model $y=a\times e^{-b\times x}+c$ is assumed also using the nls function.

3 Results and discussion

The model output is shown in Table 1. All results are derived from model simulations. To quantify the influence of consumer-level bioconcentration and bioaccumulation of Hg²⁺ on MMHg⁺ bioaccumulation, the model was run under scenarios with and without bioaccumulation of Hg²⁺ and with and without consumer-level bioconcentration of MMHg⁺. The  % bioconcentrated is calculated as $\text{Bioconcentrated (\%)}=\frac{\text{Bioconcentrated}}{\text{Bioaccumulated}}\cdot \mathrm{100}\mathit{\%}$ and the difference (%) is calculated as $\text{Difference (\%)}=\frac{\text{Scenario}}{\text{Default}}\cdot \mathrm{100}$. The italic values in the difference category indicate when the scenario causes a change larger than 10 %. The values are based on the last 10 years of the simulation and the top 20 m of the water column, to create an average value that we can compare between the setups.

3.1 Bioaccumulation of Hg²⁺

The effect of Hg²⁺ bioaccumulation on the bioaccumulation of MMHg⁺ is shown in Table 1. To quantify the influence of consumer-level bioconcentration and bioaccumulation of Hg²⁺ on MMHg⁺ bioaccumulation, the model was run under scenarios with and without bioaccumulation of Hg²⁺ and with and without consumer-level bioconcentration of MMHg⁺. The differences are low between 1 % and −6 %. This is statistically evaluated, and the results are shown in Table 2. Wilcoxon's signed rank test shows that bioaccumulation of Hg²⁺ has no significant impact on the bioaccumulation of MMHg⁺ (p=0.67). Furthermore, the Bayesian t-test shows with a BF₁₀=0.40 that the results are 2.5 times more likely under the null hypothesis of no effect than under the alternative hypothesis.

Table 1 The bioaccumulated MMHg⁺, the percentage of bioaccumulated MMHg⁺ that originates from bioconcentraton, and the bioaccumulated MMHg⁺ in the scenario without bioaccumulation of Hg²⁺ and the bioconcentration of MMHg⁺ in consumers and the difference to the default scenario. The % bioconcentrated is calculated as $\text{Bioconcentrated (\%)}=\frac{\text{Bioconcentrated}}{\text{Bioaccumulated}}\cdot \mathrm{100}$ and the difference (%) is calculated as $\text{Difference (\%)}=(\mathrm{1}-\frac{\text{Scenario}}{\text{Default}})\cdot \mathrm{100}$. Percentages above 10 % are highlighted in italic.

Table 2 The results of the statistical test performed to evaluate the difference between the scenarios and the base case. The high p-value (p=0.67) and below 1 Bayes Factor (BF₁₀=0.40) indicate that there is no significant difference between the mean of the base case and the scenario without Hg²⁺ bioaccumulation and that the change that the mean is equal is 2.5 times larger than the chance that there it is not. The difference between the scenario without MMHg⁺ bioconcentration is significant (p<0.001) and the change that there is a difference in the mean bioaccumulation of MMHg⁺ caused by the bioconcentration of MMHg⁺ in consumers is 5.96 times higher than the change that there is no difference in the mean.

3.2 Bioaccumulation of MMHg⁺

The MMHg⁺ bioaccumulation for all biota functional groups in the different setups and scenarios and the percentage of bioaccumulated MMHg⁺ originating from bioconcentration are shown in Table 1. These results show that the relative contribution of bioconcentration to the MMHg⁺ content is low in microzooplankton (4 %–6 %) while it is higher in mesozooplankton (5 %–11 %) higher in fish 1 (13 %–22 %), while lower in fish 2 (8 %–14 %) and higher in macrobenthos (14 %–27 %). The relative contribution of direct bioconcentration to the MMHg⁺ bioaccumulation in zooplankton, especially microzooplankton, is lower than in higher trophic levels of animals. In our model, this occurs because of the extremely high turnover rate of zooplankton. This “grow fast, die young” approach results in less MMHg⁺ bioconcentration with higher relative contributions due to feeding caused by the high feeding rate of zooplankton.

In longer-lived fish, we see higher contributions of bioconcentration. Although these contributions are higher, they align with the experiments of (Wang and Wong, 2003) and the observations of 15 % by Hall et al. (1997). Both fish 1 and fish 2 have the same bioconcentration and release rates, so it is in line with expectations that the relative contribution of direct bioconcentration in fish 2 is lower than in fish 1, since it gets more MMHg⁺ from its higher trophic level diet.

There is a great difference in the importance of bioconcentration of MMHg⁺ in macrobenthos between the Southern and Northern North Sea. This difference is especially notable in the direct bioconcentration in macrobenthos, which is 27 % of the total bioaccumulated MMHg⁺ in the Northern North Sea and only 14 % in the Southern North Sea. This difference is caused by the low intake of MMHg⁺ from food by macrobenthos in the Northern North Sea. Since the water column is stratified during spring and summer, macrobenthos cannot directly feed on the phyto- and zooplankton bloom. Because of this, they are dependent on sinking detritus. The detritus has a lower MMHg⁺ content than living material and consequently, the MMHg⁺ intake in Northern North Sea macrobenthos is lower, and thus the relative importance of bioconcentration is higher.

The percentage differences between the base case and the alternative setups are also shown in Table 1. Notably, the change in MMHg⁺ concentrations is substantially larger in the scenario without consumer-level bioconcentration. The greatest difference between the base case and the setup without consumer-level bioconcentration occurs in fish 2, with increases ranging from 28 % to 49 %, while this effect is smaller in lower-trophic level biota. Interestingly, in the Southern North Sea, a reduction in the bioaccumulation of MMHg⁺ in primary producers is observed in both the scenario without Hg²⁺ bioaccumulation (3 %–6 %) and the scenario without consumer-level MMHg⁺ bioconcentration (3 %–7 %). This reduction is likely linked to the removal of bioaccumulated Hg from the water column by macrobenthos feeding on pelagic food sources; consequently, it is not present in the two deep-water setups. Nevertheless, the overall effect on MMHg⁺ bioaccumulation remains limited between −2 % and 2 %.

3.3 Bioaccumulation of MMHg⁺ and trophic level

The relationship between trophic level and MMHg⁺ bioaccumulation is plotted in Fig. 1a. Since we assume biomagnification to be an exponential effect per trophic level on top of bioconcentration, the model is fitted as an exponential function with the average MMHg⁺ content of phytoplankton as the origin. Figure 1b expands on this by showing the relationship between the trophic level and the effect of both investigated drivers on the bioaccumulation of MMHg⁺. This shows that while there is no effect of the bioaccumulation of Hg⁺ on the bioaccumulation of MMHg⁺, the setup without bioconcentration in consumers has 15 % less MMHg⁺ bioaccumulation per trophic level than the setup which does include this interaction.

Figure 1 Panel (a) shows the trophic level vs MMHg⁺ bioaccumulation for the base case and the 2 scenarios. The base case and scenario “No Hg²⁺ bioaccumulation” have the relationship $\mathrm{0.0075}+\mathrm{6}\times {\mathrm{10}}^{-\mathrm{4}}\times e^{\mathrm{1.2}\cdot \text{TrophicLevel}}$ and the scenario “No MMHg⁺ bioconcentration” has $\mathrm{0.0075}+\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}\cdot e^{\mathrm{1.09}\cdot \text{TrophicLevel}}$. Panel (b) shows that there is an exponential increase in MMHg⁺ with trophic level, which is higher for the base case and the scenario without Hg bioaccumulation than for the scenario without MMHg⁺ bioconcentration. Panel (b) expands on this and demonstrates no significant effect on the importance of trophic level on the effect of the bioaccumulation of Hg²⁺ on MMHg⁺ bioaccumulation. There is a reduction of 15 % per trophic level caused by the bioconcentration of MMHg⁺.

3.4 Seasonality of the difference in MMHg⁺ bioaccumulation

The seasonality of the difference in MMHg⁺ bioaccumulation caused the bioaccumulation of Hg²⁺ and the bioconcentration of MMHg⁺ in consumers is shown in Fig. 2. For each calendar day (1 January, 2 January, etc.), the modeled daily values from each of the last 10 years of the simulations were averaged. The resulting time series represents an annual cycle of average daily conditions. From the producers' functional groups, only the diatoms are shown as the reaction is not group-specific but rather caused by changes in dissolved Hg²⁺ and MMHg⁺ which means the difference caused for all phytoplankton groups was the same. This shows that, while the scale depends on the setup, there are interactions that consistently occur. In low trophic levels, such as phytoplankton and microzooplankton, the bioaccumulation of Hg²⁺ causes a seasonal response in the MMHg⁺ bioaccumulation in phytoplankton, which is consequently observable in low trophic level biota such as microzooplankton. While this reduction in MMHg⁺ would compound into higher trophic levels, its effects in higher trophic level animals dwarf in comparison to the difference caused by incorporating the bioconcentration of MMHg⁺ in consumers, and it does not cause a difference larger than 3 % in either fish 1 or fish 2 in any of the setups.

Figure 2 The seasonality of the difference in the bioaccumulation of MMHg⁺ caused by the bioaccumulation of Hg²⁺ and the bioconcentration of MMHg in consumers for (a) the Gotland Deep, (b) the Northern North Sea and (c) the Southern North Sea. In high-trophic-level such as fish 1 and fish 2 there is low seasonality and the effect of the bioconcentration of MMHg⁺ in consumers is high while the effect of the bioaccumulation of Hg⁺ is low. In low trophic levels, notably diatoms and microzooplankton there is strong seasonal component. The bioaccumulation of MMHg⁺ is up to 5 % lower in diatoms in the Southern North Sea if the bioaccumulation of Hg²⁺ is modeled in late summer when biomass is high. But the bioaccumulation of Hg²⁺ does not lead to a notable (>5 %) difference at any moment in fish.

4 Sensitivity analyses

4.1 Sensitivity of the consumer bioconcentration rate

The results of the first sensitivity study, in which the bioconcentration rate of consumers is altered, are shown in Fig. 3. Figure 3a illustrates that the MMHg⁺ contribution from bioconcentration in consumers is linearly related to the consumer bioconcentration rate scaling factor. Thus, in bioaccumulation modeling, altering the bioconcentration rate by half or double yields the same relative effect on fish 2's MMHg⁺ content from direct bioconcentration. Based on Table 1, we can see that in the Gotland Deep, the difference between the simulation with and without consumer bioconcentration is 0.0183 ng Hg mgC⁻¹. This means that parameterizing a bioconcentration rate double the real rate would result in a 0.0183 ng Hg mgC⁻¹ overestimation of MMHg⁺ bioaccumulation in fish 2, while selecting bioconcentration rates half the true values would result in a reduction of 0.00915 ng Hg mgC⁻¹. However, the relative contribution of bioconcentration to total MMHg⁺ bioaccumulation follows a non-linear pattern, as shown in Fig. 3b. This non-linearity occurs because the total MMHg⁺ in fish 2 is influenced by both bioconcentration in consumers and bioconcentration in producers. When the consumer bioconcentration scaling factor is 0, bioconcentration in consumers makes no contribution to fish 2's MMHg⁺ levels, thus the percentage difference is also 0. However, the contribution of consumer-level bioconcentration can never reach 100 %, because bioconcentration in producers and consequent biomagnification from lower trophic levels always contributes to the total MMHg⁺ burden in fish. In the same way as in the results shown in Table 1, the relative importance of bioconcentration is consequently highest in the Northern North Sea, followed by the Southern North Sea and lowest in the Gotland Deep.

Figure 3 (a) show the effect of the bioaccumulation of MMHg⁺ per trophic level in the Gotland Deep. This shows an increase 0.0036±0.00010 ng Hg mg⁻¹ in fish 2 MMHg⁺ bioaccumulation for every 0.2 step increase in the consumer bioconcentration scacling factor. (b) shows the percentage difference due to bioaccumulation in fish 2 with different consumer bioconcentration scaling factors in all setups. GD refers to the Gotland Deep, NNS to the Northern North Sea and SNS to the Southern North Sea. When the consumers bioconcentration scaling factor is 0, the percentage difference due to bioconcentration is 0 %. As this increase the percentage increases. The relationship between the consumers bioconcentration factor and the percentage difference due to consumer bioconcentration is plotted assuming an saturating exponential relationship.

4.2 Sensitivity of the producer bioconcentration rate

The results of the second sensitivity study are shown in Fig. 4. Here, rather than the consumer bioconcentration rate, the producers' bioconcentration rates are multiplied by a scaling factor. Again, the effect of this scaling on the bioaccumulation in all trophic levels is visualised in Fig. 4a, and the effect of this scaling on the relative importance of consumer bioconcentration on MMHg⁺ bioaccumulation is shown in Fig. 4b. If the bioconcentration scaling factor is 0, there is still MMHg⁺ bioaccumulation in fish 2, both from direct bioconcentration and from bioconcentration in consumers and consequent biomagnification. The increase in fish 2 MMHg⁺ per step of 0.2 in the scaling factor is 0.0083±0.00030 ng Hg mg⁻¹. The relative contribution of consumer bioconcentration on MMHg⁺ bioaccumulation is shown in Fig. 4b. An important note here is that while we scaled the bioconcentration factor of producers and consumers, MMHg⁺ can also be bioaccumulated via the partitioning to DOM detritus and consequent biomagnification. This is especially important in the Northern North Sea. In the seasonally stratified water column, macrobenthos cannot feed directly off the phytoplankton bloom; thus, the dying and sinking of particles is an important flux that is consumed by the benthos. Benthos, in turn, is an important food source for fish 2. So scaling the producer bioconcentration rate has less effect in the Northern North Sea. In the Gotland Deep, the opposite is true; because the deep water is anoxic, there is no macrobenthos in the model. This means that the entire ecosystem is pelagic and detritus is less important than direct consumption of the phytoplankton bloom.

Figure 4 (a) illustrates the influence of scaling the producers bioconcentration rate of MMHg⁺ on the MMHg⁺ bioaccumulation at each trophic level in the Gotland Deep. This shows an increase of 0.0084±0.00032 ng Hg mg⁻¹ in fish 2 MMHg⁺ with every 0.2 increase in the producers scaling factor. (b) shows importance of consumer-level bioconcentration across all setups by showing the percentage difference in the bioaccumulation of MMHg⁺ in fish 2 between the setup with and without consumer-level bioconcentration. This shows that in all cases the percentage difference is high when the producer bioconcentration factor is low, and that this percentages decreases with an increasing producer bioconcentration scaling factor.

5 Evaluation of the hypotheses

Our results show that the MMHg⁺ content of high-trophic-level fish is a combination of the direct uptake of biomagnification by these high-trophic-level fish and the MMHg⁺ they take up via their diet. The MMHg⁺ content of their diet is, in turn, made up of the bioconcentration in each trophic level, including producers and consumers.

5.1 Evaluation hypotheses 1; the effect of Hg²⁺ bioaccumulation on MMHg⁺ bioaccumulation

Based on the results of the statistical analysis shown in Table 2, we can see that there is no significant difference (p=0.67) caused by Hg²⁺ bioaccumulation on MMHg⁺ bioaccumulation. Based on the Bayesian t-test, we estimate that the change is $\mathrm{1}/\mathrm{0.40}=\mathrm{2.5}$ times greater than the mean, which is equal to that it is not. We do note that the seasonal changes in the total Hg concentration due to bioaccumulation change the bioaccumulation of MMHg⁺ at the base of the food web, and we can see this change in phytoplankton and low trophic level consumers, but it does not cause a notable (>5 %) change in MMHg⁺ bioaccumulation in fish. Based on these results, we conclude that Hg²⁺ bioaccumulation does not play a major direct role in the bioaccumulation of MMHg⁺ in our model.

5.2 Evaluation hypotheses 2; the effect of MMHg⁺ bioconcentration in consumers on MMHg⁺ bioaccumulation

Based on the statistical results shown in Table 2 we conclude that there is a significant difference between the base case and the scenario without consumer bioconcentration (p<0.001). We base our conclusion on the Bayesian t-test that the chance that the mean bioaccumulation is different without consumer-level bioconcentration is 5.96 times larger than the chance that there is no difference in the mean. Based on the results, we conclude that the bioconcentration of MMHg⁺ in consumers is a significant contributor to the bioaccumulation of MMHg⁺.

6 Model limitations

6.1 Scope of the biogeochemical model

There are some limitations to our model. First, estimates of the biomagnification factor of MMHg⁺ range between 2–10. Our model represents the estimations of the lower end. The bioconcentration factor is probably more important in low biomagnification food webs. Another limitation is that our model stops at trophic level 3.7. This is a high-trophic level that can represent piscivorous species, but many marine species that are consumed by humans, such as tuna, great marlin, and cod, can have higher trophic levels of up to 4.8 (Nilsen et al., 2008; Sarà and Sarà, 2007).

When the absolute concentration of MMHg⁺ increases at higher trophic levels, the relative increase in the importance of direct bioconcentration per trophic level likely decreases. Our modeled top predator with a trophic level of 3.7 has a high-trophic-level for a coastal ocean, but there is a marine biota with even higher trophic levels in our model domain, such as marine mammals. Without a dedicated modeling study to simulate the diet and bioconcentration of even higher trophic levels, we cannot simply extrapolate our findings to predict the importance of MMHg⁺ bioconcentration in their MMHg⁺ bioaccumulation.

6.2 Uncertainty in parameters

Overall, the most important driver of our model is the fraction of MMHg⁺ that is bioaccumulated by bioconcentration for each trophic level, as this drives the relative importance of bioconcentration at the higher trophic levels. This is influenced by both the direct bioconcentration in that trophic level and the bioconcentration of animals in lower trophic levels. This relationship is quantified in the sensitivity study. The contribution of bioconcentration in zooplankton of 4 %–11 % is in line with the < 20 % reported by Schartup et al. (2018), and the contribution of bioconcentration in fish between 8 %–22 % is in line with the study by Wang and Wong (2003).

The main uncertainty for the fraction of MMHg⁺ that originates from bioconcentration is the parameterization of bioconcentration and biomagnification. Both the bioconcentration and biomagnification of zooplankton are based on the work of Tsui and Wang (2004) on water fleas (Daphia Pulex) and for fish this is based on Wang and Wong (2003) and their work on the Indo-Pacific species sweetlips. Although water fleas are common in the Baltic Sea, they are not in the North Sea, and sweetlips live neither in the North nor in the Baltic Sea. This means that the most important parameters in our model are not based on the animals they represent in our model. Although it is unfeasible to have dedicated bioaccumulation studies in every animal or functional group, there are currently not enough studies to verify whether these rates would differ between the circumstances in our model and those in the experiment. Drivers that might influence these factors are the size of the biota, physical circumstances such as temperature and salinity, or if there is a seasonal effect related to the activity of the animals. It would greatly improve our ability to model the bioaccumulation of MMHg⁺ if more information on these different drivers were available.

6.3 Limitations due to model design

It is a deliberate choice to perform this study in 1D idealized water column models, as it allows a clear overview of the driving processes and enables us to generalize our findings. In this way, we can provide a general conclusion based on the biomagnification and bioconcentration rates of the biota that are presented in laboratory and field studies. However, it poses limitations compared to real fish by omitting spatial variability. Locally variable circumstances, such as the seasonally dependent flow of Hg from rivers to the ocean, could cause the importance of bioconcentration on MMHg⁺ bioaccumulation to be regionally different.

Although the model can predict the importance of MMHg⁺ bioconcentration, it cannot evaluate the importance of the bioconcentration of gaseous Hg species, such as Hg⁰ and DMHg. These gaseous Hg species are assumed not to biomagnify because they are not polar but could bioconcentrate. Because the gills of fish are optimized to facilitate the exchange of gases between water and fish blood, these gaseous Hg species can likely bioconcentrate in organisms, but it is unclear to what degree. To model and evaluate the importance of this interaction, studies must be performed first, investigating both the bioconcentration and release rates of these gaseous species and their fate in the organism. In particular, the effect of the bioconcentration of DMHg on the concentration of MeHg at higher trophic levels will be very dependent on whether DMHg stays gaseous in the organism and is excreted quickly via the gills, or whether DMHg is demethylated in MMHg⁺ and further biomagnified in the food chain. Since DMHg concentrations are low in the North and Baltic seas, this is unlikely to play a major role in our setups, but it could influence the importance of bioconcentration on the MMHg⁺ content of higher trophic level fish in seas with higher DMHg concentrations, such as the open oceans and the Mediterranean Sea.

Figure 5 While the direct contribution of bioconcentration in consumers is limited, its cumulative effect increases with trophic level. The background and organism illustrations were generated using GPT-4.1 and manually combined, with labels added afterwards.

6.4 Uncertainty of the conclusion

The results of our model represent just one possible outcome based on a regional setup representing the North and Baltic Seas, and the importance of bioconcentration can vary greatly depending on the bioconcentration factors of all species in the trophic chain. We can assess the expected range of importance of consumer-level bioconcentration by developing theoretical maximum and minimum values based on observational studies. We can estimate that direct bioconcentration in zooplankton may account for up to 50 %, based on Lee and Fisher (2017), and similarly for mid-trophic level fish, based on Wang and Wong (2003).

We can use this to estimate the maximum expected contribution of consumer-level bioconcentration on bioaccumulation by making two assumptions: (1) bioconcentration in both copepods and fish lies between 0 % and 50 % and is equal across all trophic levels, and (2) the food chain is linear, meaning that trophic level 3 feeds exclusively on trophic level 2, which feeds exclusively on trophic level 1. Under these assumptions, we can estimate the percentage of MMHg⁺ in the diet of a given trophic level that originated from bioconcentration in primary producers as:

$$\begin{array}{}\text{(17)}& {\text{PBC\hspace{0.17em}\%}}_n=(\mathrm{1}-\text{BC}{)}^{n-\mathrm{1}}\times \mathrm{100}\mathit{\%}\end{array}$$

where:

• PBC %_n is the percentage of MMHg⁺ in the diet of trophic level n that originates from bioconcentration at the primary producer level,

• BC is the fraction (0–1) of MMHg⁺ at each trophic level originating from bioconcentration.

Although this is a simplification, it illustrates that a high bioconcentration estimate of 50 % results in only 12.5 % of MMHg⁺ in the diet of a trophic level 4 fish originating from bioconcentration in primary producers, meaning that 87.5 % originates from consumer-level processes. Even a low estimate of 10 % still results in 27.1 % of MMHg⁺ in the diet of the same high-trophic-level fish originating from consumer-level bioconcentration.

The degree to which this interaction contributes to overall bioaccumulation depends on numerous additional factors that are not yet fully understood, including the size distribution of phytoplankton at the base of the food web, the trophic structure, consumer metabolic and respiration rates, and the assimilation efficiency of MMHg⁺ from the diet. This complexity makes it difficult, if not impossible, to provide a definitive estimate of the importance of consumer-level bioconcentration and the uncertainty of the interaction. However, based on the bioconcentration rates provided in the current literature, we conclude that this process plays a key role in the bioaccumulation of MMHg⁺ in higher trophic levels.

7 Conclusions

Our results show that the bioaccumulation of MMHg⁺ in our model with and without the bioaccumulation of Hg²⁺ is not significantly different, while this is the case for the model with and without the bioconcentration of MMHg⁺. We show that the bioconcentration of MMHg⁺ in consumers becomes more important at higher trophic levels because it is an effect of the sum of all trophic levels before it. Our main conclusion is summarized in Fig. 5 and shows that while the direct bioconcentration direct bioconcentration only accounts for 8 %–14 % of MMHg⁺ bioaccumulation in our highest trophic level fish, the total effect of bioconcentration in consumers accounts for 28 %–49 %. This effect increases with the trophic level and the percentile contribution of the cumulative effect of MMHg⁺ bioconcentration in consumers on MMHg⁺ bioaccumulation is 15 % per trophic level.

Because of this, we reject the first hypothesis that bioaccumulation of Hg²⁺ lowers MMHg⁺ bioaccumulation and accept our second hypothesis that bioconcentration of MMHg⁺ increases bioaccumulation of MMHg⁺ in higher trophic level fish. We supplement the second hypothesis by quantifying the effect as an average increase in bioaccumulated MMHg⁺ of 15 % per trophic level.

These results demonstrate that to model the bioaccumulation of MMHg⁺, the bioaccumulation of Hg²⁺ can be ignored to save computational resources. However, the bioconcentration of MMHg⁺ on the other hand, is an important interaction that should be taken into account. When modeling MMHg⁺ bioaccumulation at higher trophic levels, the inclusion of MMHg⁺ bioconcentration is therefore essential.