Articles | Volume 19, issue 15
Biogeosciences, 19, 3663–3682, 2022
Biogeosciences, 19, 3663–3682, 2022
Research article
05 Aug 2022
Research article | 05 Aug 2022

Assessing the spatial and temporal variability of methylmercury biogeochemistry and bioaccumulation in the Mediterranean Sea with a coupled 3D model

Assessing the spatial and temporal variability of methylmercury biogeochemistry and bioaccumulation in the Mediterranean Sea with a coupled 3D model
Ginevra Rosati, Donata Canu, Paolo Lazzari, and Cosimo Solidoro Ginevra Rosati et al.
  • National Institute of Oceanography and Applied Geophysics – OGS, Trieste, 34010, Italy

Correspondence: Ginevra Rosati (


Human exposure to mercury (Hg) is a cause of concern, due to the biomagnification of the neurotoxic species monomethylmercury (MMHg) in marine ecosystems. Previous research revealed that commercial fish species in the Mediterranean Sea ecosystems are particularly enriched in Hg, due to a combination of physical and ecological factors. Since the fate of Hg depends on the interactions among several biogeochemical and physical drivers, biogeochemical modeling is crucial to support the integration and interpretation of field data. Here, we develop and apply a coupled transport–biogeochemical–metal bioaccumulation numerical model (OGSTM–BFM–Hg) to simulate the biogeochemical cycling of the main Hg species (HgII, Hg0, MMHg, and DMHg) in seawater, organic detritus, and through the planktonic food web. The model is applied to a 3D domain of the Mediterranean Sea to investigate the spatial and temporal variability of methylmercury (MeHg) distribution and bioaccumulation and major uncertainties in Hg cycling. Model results reproduce the strong vertical and zonal gradients of MeHg concentrations related to primary production consistently with the observations and highlight the role of winter deep convection and summer water stratification in shaping the MeHg vertical distribution, including subsurface MeHg maximum. The modeled bioaccumulation dynamics in plankton food webs are characterized by a high spatial and temporal variability that is driven by plankton phenology and is consistent with available field data of Hg concentrations in plankton, as well as with other indicators, such as bioconcentration factors (BCFs) and trophic magnification factors (TMFs). Model results pointed out that the increment in water temperature linked to a decline of deep convection can cause an increase in water MeHg concentrations with cascading effects on plankton exposure and bioaccumulation.

1 Introduction

The anthropogenic alteration of the mercury (Hg) biogeochemical cycle has led to global enrichment of Hg concentrations in all the environmental compartments (Amos et al., 2015; UNEP, 2019; Zhang et al., 2014a, b) and concerns over human exposure to neurotoxic monomethylmercury (MMHg), which is produced in the ocean and biomagnifies in marine food webs (UNEP, 2019). The ocean has absorbed 50 % of anthropogenic Hg historical emissions, 35 % of which are currently stored in the water and the rest in the sediments (Zhang et al., 2014b): relative to pre-anthropogenic levels, the Hg enrichment is 230 % in surface ocean waters, 25 % in intermediate waters, and 12 % in deep waters (UNEP, 2019).

Models and observations (Lamborg et al., 2014; Zhang et al., 2014b) suggest that the oceanic water column contains 280–370 Mmol of anthropogenic Hg (equivalent to about 70 000 t), a large part of which is recycled within surface and intermediate water due to the biological carbon pump and the microbial loop (e.g., Zhang et al., 2018). Hg is scavenged from surface water by organic particles through adsorption and uptake, transported downward via particles sinking, and released in the dissolved phase following particle degradation. The locally increased microbial activity and availability of inorganic Hg at the depth of particle remineralization are also thought to promote Hg methylation (Cossa et al., 2009; Heimbürger et al., 2010; Sunderland et al., 2009). Concurrent biological and photochemical transformations occur in the water column, driving the interconversions among the main mercury species in seawater (inorganic oxidized Hg, HgII; monomethylmercury, MMHg; elemental Hg, Hg0; and dimethylmercury, DMHg). Other important processes driving a fast cycling of Hg species in surface water are the atmospheric deposition of oxidized Hg and the gaseous exchange of Hg0 (Jiskra et al., 2021). The resulting vertical distribution of Hg species mirrors the local equilibrium among a number of biogeochemical and physical drivers. However, recent comprehensive reviews of the current understanding of Hg fate in the ocean (Bowman et al., 2020; Mason et al., 2012; Sonke et al., 2013) and in the Mediterranean Sea (Cossa et al., 2022) highlighted still high uncertainty in describing and quantifying most of these processes.

Methylmercury (MeHg; defined as the sum of MMHg and DMHg) is more frequently measured than the individual methylated species. In oceanic and Mediterranean areas with high primary production, the MeHg vertical distribution is characterized by low concentrations at the surface due to photodegradation and concentration maxima coincident with the maximum in apparent oxygen utilization, a proxy for heterotrophic activity (Cossa et al., 2009; Sunderland et al., 2009). These observations suggest that most oceanic MeHg is produced in situ. This conclusion is also supported by estimates of MeHg fluxes in the ocean (Mason et al., 2012; Zhang et al., 2020) and by studies based on stable Hg isotopes (Blum et al., 2013; Motta et al., 2019). The role of DMHg in controlling MMHg dynamics remains puzzling: recent observations (Cossa et al., 2017) from the Northwestern Mediterranean Sea support the idea that DMHg is produced from HgII in highly productive waters and then degraded to MMHg (Conaway et al., 2009). However, only an abiotic pathway has been directly identified for the formation of DMHg (Jonsson et al., 2016), and rates of DMHg formation from HgII measured in polar water were negligible compared to methylation rates from HgII to MMHg (Lehnherr et al., 2011).

Tunas and other commercial fish species from the Mediterranean Sea are particularly enriched in Hg in comparison to the same species from other geographical areas (Cossa et al., 2012; Cossa and Coquery, 2005; Harmelin-Vivien et al., 2009; Tseng et al., 2021). It has been suggested that, along with other ecological features of the whole food web, this is due to a shallower occurrence of the MeHg concentration maxima compared to the ocean which results in higher phytoplankton exposure and bioaccumulation (Cossa et al., 2012). Available observations on MeHg distribution in the Mediterranean Sea suggest that concentrations are lower in the oligotrophic waters of the Ionian Sea than in the mesotrophic waters of the Northwestern Mediterranean (Cossa et al., 2009, 2022) and Adriatic Sea (Kotnik et al., 2015).

Recent efforts to couple the dynamics of biological carbon pump and microbial loop with Hg dynamics, including bioaccumulation in the lower food web, provided global-scale assessments of MeHg production and bioaccumulation (Wu et al., 2020, 2021; Zhang et al., 2020). The only study modeling the Hg dynamics in the Mediterranean Sea (Žagar et al., 2007) did not include any biological component. Much of the current knowledge on the dynamics of MeHg in the open sea was acquired only afterwards (Bowman et al., 2016, 2015; Cossa et al., 2017, 2009; Heimbürger et al., 2010; Lehnherr et al., 2011; Monperrus et al., 2007; Munson et al., 2015, 2018; Ortiz et al., 2015; Sunderland et al., 2009).

Here, we couple a model for Hg biogeochemistry (Canu et al., 2019; Rosati et al., 2018, 2020; Zhang et al., 2020) with the 3D transport biogeochemical model OGSTM–BFM (Cossarini et al., 2021; Lazzari et al., 2010, 2021; Salon et al., 2019) to investigate spatial and temporal variability of MeHg dynamics in the Mediterranean Sea. The coupled model, which has a 1/16 horizontal resolution, describes the biotic and abiotic interconversions among four Hg species (HgII, MMHg, Hg0, and DMHg), as well as other key processes of the Hg cycle such as partitioning to organic detritus, sinking to the seabed, exchange of gaseous Hg with the atmosphere, and bioaccumulation in four phytoplankton and four zooplankton functional groups (Lazzari et al., 2012; Vichi et al., 2015). The agreement with observations has been improved by tuning the sinking velocity of organic detritus and by exploring the sensitivity to alternative parameterizations of Hg methylation. A first assessment of the uncertainties concerning the fate of riverine Hg loads is also achieved and discussed along with other terms of the Hg budget in the Mediterranean Sea. Model results are compared to available observations discussing the implications of spatial and temporal variability of modeled MeHg distributions in water and plankton.

2 Methods

2.1 Study area

The Mediterranean Sea (Fig. 1) is a semi-enclosed basin characterized by decreasing west–east gradient of primary production resulting from the superposition of biological pump and inverse estuarine circulation, as well as from the limited impact of riverine loads on the open-ocean dynamics (Crise et al., 1999, Crispi et al., 2001). Due to the formation of intermediate and deep water generating thermohaline circulation, the Mediterranean is considered a miniature of the ocean (Pinardi et al., 2019). The large-scale basin circulation (Pinardi et al., 2015, 2019; Pinardi and Masetti, 2000, and references therein) is composed of three thermohaline cells (Fig. 1). A zonal vertical circulation belt is associated with the inflow of shallow Atlantic Water (AW) at the Strait of Gibraltar that becomes progressively saltier due to evaporation moving eastward (MAW, Modified Atlantic Water) and eventually sinks at intermediate depths forming the Levantine Intermediate Water (LIW), outflowing at Gibraltar. The other cells are driven by dense water formation due to winter cooling, salting, and sinking of surface water. Western Mediterranean Deep Water (WMDW) originates in the Northwestern Mediterranean Sea, while Eastern Mediterranean Deep Water (EMDW) is formed in the Adriatic Sea but also in the Aegean and the Levantine Sea. The mixing between these deep water masses is inhibited by the shallowness of the Strait of Sicily (about 500 m). Long-term observations of ocean color (Bosc et al., 2004; Bricaud et al., 2002; D'Ortenzio and D'Alcalà, 2009) and coupled transport–biogeochemical models (Di Biagio et al., 2019, 2021, Lazzari et al., 2012, 2014, 2021) supported the characterization of different biogeoprovinces based on the seasonal, inter-annual, and high-frequency variability of biogeochemical processes in different areas (subbasins). Biogeochemical modeling also revealed the importance of subsurface plankton blooms, which are not captured by ocean color satellite observations, highlighting the existence of a deep chlorophyll maximum that becomes deeper with increasing water oligotrophy. The most productive areas are the mesotrophic western subbasins of the Alboran Sea (Alb), Northwestern and Southwestern Mediterranean (Nwm and Swm). The eastern Mediterranean is ultra-oligotrophic, apart from the Northern Adriatic Sea (Nad), which is a coastal sea, and other “intermittently blooming” areas (i.e., with high inter-annual variability) such as the Southern Adriatic Sea (Sad), the Aegean Sea (Aeg), and the Rhodes Gyre in the Levantine Sea (Lev).

Figure 1Bathymetric map of the Mediterranean Sea (from Ocean Data View, Schlitzer, 2014). The white contours indicate the subdivision of the model domain in subbasins: Alboran Sea (Alb), Southwestern Mediterranean (Swm1 and Swm2), Northwestern Mediterranean (Nwm), Tyrrhenian Sea (Tyr1 and Try2), Northern Adriatic Sea (Nad), Southern Adriatic Sea (Sad), Ionian Sea (Ion1, Ion2, and Ion3), Aegean Sea (Aeg), and Levantine Sea (Lev1, Lev2, Lev3, and Lev4). The light blue line shows the zonal cell driven by the surface inflow of Atlantic Water (AW) that moves eastward, forming Modified Atlantic Water (MAW) and Levantine Intermediate Water (LIW), outflowing at intermediate depths at the Strait of Gibraltar. The dark blue lines show the path of meridional cells related to winter convection in the Nwm, Nad, Sad, Aeg, and Lev originating the Western and Eastern Mediterranean Deep Water (WMDW and EMDW).

2.2 The coupled model

The Biogeochemical Flux Model (BFM; Vichi et al., 2015) is a model that simulates the cycling of carbon, nitrogen phosphorus, and silica through dissolved, living, and non-living particulate phases of the marine environment. The description of the planktonic food web is based on the plankton functional type (PFT) approach: the pool of species having similar traits (nutrients affinities, light affinity, prey–predator relationships) are described by a single variable. In the BFM, PFTs have variable stoichiometry and simulate the lower trophic levels of marine food webs up to carnivorous zooplankton (Sect. 2.2.1). The OGSTM–BFM model has been coupled offline to the ocean general circulation model (NEMO) to investigate different issues related to the Mediterranean biogeochemistry (Canu et al., 2015; Cossarini et al., 2015, 2021; Lazzari et al., 2021, 2012, 2014, 2016), and it is now routinely used as the backbone of the biogeochemical component of the Copernicus CMEMS Mediterranean system (, last access: 25 July 2022, Salon et al., 2019; Terzić et al., 2019). The integration of a Hg cycling module in the OGSTM–BFM model adds more flexibility to this system and opens the way to novel operational oceanography applications. Here the OGSTM–BFM–Hg model is implemented (Fig. 2) by coupling the OGSTM–BFM model with the equations describing mercury (Hg) species partitioning, transformation, transport (Sect. 2.2.2), bioconcentration in phytoplankton (Sect. 2.2.3), and trophic transfer to zooplankton (Sect. 2.2.4). The coupled model is used to investigate the Hg cycling and bioaccumulation dynamics (Sect. 2.2.5) at the base of the Mediterranean Sea food webs. The model has a 1/16 horizontal resolution (approximately 6 km) and 72 vertical levels and runs with a 300 s time step. A model spin-up of 13 years forced with climatological data (Sect. 2.2.6) is performed, and the results of the last year of simulation (2017) are analyzed.

Figure 2Overview of the coupling in the OGSTM–BFM–Hg model, which integrates the biogeochemistry and bioaccumulation of Hg species in the OGSTM–BFM model, previously coupled to the transport model NEMO and to a multispectral biotical model. The thick grey arrows with text highlight variables and fluxes read by the Hg biogeochemical model (Fig. 3).

2.2.1 BFM model dynamics

The BFM model dynamically simulates nine plankton functional types (PFTs) representative of four phytoplankton groups (picophytoplankton, nanophytoflagellates, diatoms, and large phytoplankton), four zooplankton groups (heterotrophic nanoflagellates, microzooplankton, omnivorous mesozooplankton, and carnivorous mesozooplankton), and one group of heterotrophic bacteria (Vichi et al. 2015). The model reproduces the fluxes of carbon, nitrogen, phosphorous, and silicates from the inorganic nutrient pool to organisms (phytoplankton, zooplankton, and bacteria) and organic compounds pools (particulate and dissolved organic carbon, POC and DOC) with a variable stoichiometric formulation. Non-living organic matter includes one class of organic detritus (POC) and three classes of DOC (labile, semi-labile, and semi-refractory).

2.2.2 Biogeochemical dynamics in the Hg model

The Hg biogeochemical model (Canu and Rosati, 2017; Canu et al., 2019; Melaku Canu et al., 2015; Rosati et al., 2018, 2020, 2022b, Zhang et al., 2020) simulates the cycling of four mercury species in marine water (Fig. 3): inorganic divalent Hg (HgII), monomethylmercury (MMHg), elemental mercury (Hg0), and dimethylmercury (DMHg). Equations (1)–(4) describe the time variation of non-conservative tracers representing Hg species, omitting advective and diffusive processes that are resolved by the transport model.


Field observations suggest that particulate Hg dynamics are controlled by regenerative scavenging mechanisms dependent on adsorption, remineralization, and sinking, with little influence of desorption (Lamborg et al., 2016). In the coupled OGSTM–BFM–Hg model, the dynamics POC and DOC including sinking and degradation are spatially resolved by the BFM model, and their concentrations are read by the Hg model (Figs. 2 and 3) to calculate the partitioning of HgII and MMHg between particulate species (associated with organic detritus, i.e., HgPOCII and MMHgPOC) and dissolved species (associated with dissolved organic carbon, i.e., HgDOCII and MMHgDOC, and ionic species, i.e., HgCln and MMHgCl). By combining POC and DOC concentrations with chemical partition constant KD (Table S1) (Choe et al., 2003; Choe and Gill, 2003; Lamborg et al., 2016), the fraction of each Hg species with respect to its total (fHg-POC,fHg-DOC,fHgCl,fMeHg-POC,fMeHg-DOC,fMegCl) is estimated (e.g., Eqs. 5–7 for HgII) assuming repartition at thermodynamic equilibrium. The fractions are used to calculate concentrations of particulate and dissolved Hg species from the total pool of HgII and MMHg (e.g., Eq. 8). The pool of dissolved gaseous mercury species (DGM, i.e., Hg0 and DMHg) is assumed to be entirely in the dissolved phase.


Once the partitioning is determined, the model computes the transformation fluxes based on first-order rates (kx, d−1 Table S2), the bioconcentration (Sect. 2.3) and bioaccumulation fluxes (Sect. 2.4), the sinking of HgPOC and MMHgPOC, the advective transport within the model domain and at the model boundaries (i.e., Strait of Gibraltar), and the exchange of Hg0 with the atmosphere. Microbial Hg methylation (kmet HgII), reduction (kbioredHgDOCII+HgCln), and oxidation (kbiox Hg0) are parameterized as functions of particulate organic matter remineralization (Zhang et al., 2020, 2014), which is read from the BFM model. The direct conversion of HgII to DMHg is assumed to be negligible, due to the low rates observed in the only dataset available for this process (Lehnherr et al., 2011). Although several observations linked Hg methylation to oxygen consumption and organic matter remineralization in the water column, a full mechanistic understanding of methylation and demethylation is still lacking (An et al., 2019; Wang et al., 2020). Here we assumed that both reactions involve the entire HgII and MMHg pools (Ortiz et al., 2015; Schaefer and Morel, 2009; Zhang et al., 2019). Photochemical transformations, i.e., photoreduction (kphredHgDOCII+HgCln), photooxidation (kphox Hg0), and photodegradation of MMHg (kphdem {MMHgDOC +MMHgCl}) and DMHg (kphdm2 DMHg) are assumed to act only on the dissolved fractions and are a function of the attenuated shortwave flux (Zhang et al., 2014). MMHg photodegradation is assumed to yield both Hg0 and HgII (kphdm {MMHgDOC+MeHgCl}) due to the contrasting results in experimental studies (Luo et al., 2020).

Figure 3Hg dynamics in the coupled model OGSTM–BFM–Hg. The Hg model simulates the cycling of inorganic (HgII and Hg0) and methylated (MMHg and DMHg) Hg species resolving the partitioning of HgII and MMHg between particulate (HgPOC, MMHgPOC) and dissolved species (HgDiss, as HgDOC+ HgCln, and MMHgDiss, as MMHgDOC+ MMHgCl) based on partition coefficients (KD) and concentrations of POC and DOC, the photochemical and biological transformations (dotted colored arrows) as first-order kinetics depending on the rate constants kx, the evasion of Hg0 to the atmosphere controlled by wind speed through the volatilization rate constant kvol, and the MMHg bioaccumulation (pink arrows) driven by phytoplankton uptake in four different PFTs and trophic transfer to four zooplankton PFTs.


The gas-exchange flux of Hg0 with the atmosphere (kvolHgw0-Hgatm0/KH) uses Nightingale's parameterization (Nightingale et al., 2000) for the rate constant kvol (Eq. 9), which depends on wind speed (u10, m s−1) and on the ratio of Hg Schmidt number to CO2 Schmidt number (ScHgw/ScCO2w).

(9) k vol = 0.222 u 10 2 + 0.333 u 10 Sc Hg w Sc CO 2 w - 0.5

Current knowledge on dynamics and concentrations of DMHg in the ocean and in the atmosphere does not allow a detailed modeling representation (Baya et al., 2015; Coale et al., 2018; De Simone et al., 2014; Melaku Canu et al., 2015; Nerentorp Mastromonaco et al., 2017; Rosati et al., 2018; Zhang et al., 2020), and simulated variations of DMHg (Eq. 4) depend on photodegradation, production upon MMHg methylation, and advective transport. MMHg is taken up from the water pool by phytoplankton (MMHgphy) (Sect. 2.3), and it is released from the biotic pool back to the water column by zooplankton excretion (ExzooMMHgzooCzoo) and death (DzoozooMMHgCzoo) (Sect. 2.4). Loadings and initial and boundary conditions are presented in Sect. 2.6.

2.2.3 MMHg uptake by phytoplankton in the Hg model

The bioaccumulation model simulates the MMHg uptake by the phytoplankton PFTs of the BFM model assuming that the uptake flux (Eq. 10, nmol(Hg) m(w)-3 d−1) decreases at high DOC concentrations and increases with phytoplankton surface area : volume ratios (RSV_PFT), following the empirical relation estimated by Schartup et al. (2018) from multi-species phytoplankton uptake data (Lee et al., 2017; Lee and Fisher, 2016). The bioaccumulation model also considers the phytoplankton density (ρ= 10−12 g(w.w.)µm(w.w.)-3), wet weight biomass (Bw.w.,PFT, in g(w.w.) m(w)-3), and the water MMHg concentrations ({MMHg}w, pM). Parameter values are given in Table S3.

(10) d MMHg phy , PFT d t = 0.118 R SV _ PFT exp ( - 0.08 DOC ) 1 ρ B w . w . , PFT MMHg w 10 3

Bw.w.,PFT values were derived from the plankton biomasses expressed on a intracellular carbon content basis (CPFT, in mg(C) m(w)-3) by adopting carbon to wet weight conversion factors (Table S4) that we estimated from a dataset of experimental data for various phytoplankton and zooplankton species (Jørgensen et al., 1979; Mahlmann et al., 2008).

2.2.4 Zooplankton bioaccumulation in the Hg model

According to previous studies (Schartup et al., 2018; Zhang et al., 2020), MMHg bioaccumulation in zooplankton is mostly driven by grazing (> 80 %); therefore, the trophic transfer of MMHg from phytoplankton to zooplankton PFTs (heterotrophic nanoflagellates, microzooplankton, omnivorous mesozooplankton, and carnivorous mesozooplankton) in the model is set proportionally to the carbon (C) transfer, neglecting direct uptake from seawater. Thus, the bioaccumulation of MMHg in each zooplankton PFT (Eq. 11) depends on carbon grazing fluxes (Gzoo, PFT, mg(C) m−3 d−1) and the MMHg : C ratio of the preys that are grazed (MMHgphy,PFTCphy,PFT, nmol(Hg) mg(C)-1). MMHg release to the water pool occurs following carbon excretion (Exzoo, PFT) and zooplankton death (Dzoo, PFT, mg(C) m−3) and depends on the variable MMHg : C quota (MMHgzoo,PFTCzoo,PFT, nmol(Hg) mg(C)-1) in each zooplankton group.

(11) d MMHg zoo , PFT d t = G zoo , PFT MMHg phy , PFT C phy , PFT - Ex zoo , PFT MMHg zoo , PFT C zoo , PFT - D zoo , PFT MMHg zoo , PFT C zoo , PFT

2.2.5 Indicators for MMHg bioaccumulation and biomagnification

To compare the model content of MMHg in phytoplankton and zooplankton against experimental data, the output in nmol(Hg)m(w)-3 was converted into ng gw.w.-1 by computing the quota of MMHg (QMMHg, PFT, Eq. 12) to plankton wet weight (Bw.w.,PFT, gw.w. m(w)-3) and correcting for the molar mass of Hg (200.59 g mol−1).

(12) Q MMHg , PFT = MMHg phy , PFT B w . w . , PFT 200.59

The bioconcentration factor (BCF, L kg−1) is used to synthesize the enrichment of Hg in plankton with respect to water concentrations (Harding et al., 2018; McGeer et al., 2003), and due to the normalization for MMHg water levels, it allows better intercomparability across ecosystems than Hg concentrations. We calculated BCFs (Eq. 13) for each plankton PFT as the ratio between the quota QMMHg, PFT and water MMHg concentrations, using model output at monthly resolution spatially averaged. The logarithmic scale is often used to handle BCF values (e.g., McGeer et al., 2003).

(13) BCF = 10 6 Q MMHg , PFT MMHg w 200.59

The trophic magnification factor (TMF) is used to explore biomagnification across ecosystems by expressing the increase in Hg along with trophic levels (Alava et al., 2018; Harding et al., 2018). TMF is usually calculated in field studies of the whole food web assessing both concentrations of pollutants and trophic levels through nitrogen isotopes. Here, in agreement with other modeling studies (Wu et al., 2021, 2020; Zhang et al., 2020), TMFs were estimated from the ratio of the MMHg quota in zooplankton PFTs to that of grazed phytoplankton PFTs (Eq. 14). The analysis was carried out considering the model plankton food web as composed of a lower trophic level “herbivorous” part and an upper “omnivorous” part.

(14) TMF = Q MMHg, zoo-PFTs Q MMHg,  phy-PFTs

2.2.6 Initial conditions, boundary conditions, and forcings

Model initial conditions, boundary conditions, and physical forcings at monthly resolution are set according to Lazzari et al. (2021). Initial conditions for Hg species (Table S5) in Mediterranean subbasins (Fig. 1) are extrapolated from published data from different investigations in the Mediterranean Sea (Cossa et al., 2017; Cossa and Coquery, 2005; Ferrara et al., 2003; Horvat et al., 2003; Kotnik et al., 2007, 2015), and boundary conditions are set in agreement with the observations from the GEOTRACES GA03 meridional transect (Bowman et al., 2015). River Hg load is estimated based on dissolved Hg concentration, assuming 3 pM and 3.5 % of MMHg (Cossa et al., 2017), using a higher concentration for three rivers that are known to be highly impacted by mining (Isonzo-Soča river, 25 pM) (Hines et al., 2000), industrial activities (Po river, 6 pM) (Vignati et al., 2008), or both of them (Tiber river, 6 pM) (Lanzillotta et al., 2002; Rimondi et al., 2019), assuming MMHg is the 1 % (Hines et al., 2000). Further, a sensitivity simulation accounting for the load of particulate Hg was carried out (Sect. 2.3.3) to investigate model uncertainty regarding this source. The atmospheric Hg0 concentration is set to a constant value Hgatm0=1.6 ng m−3 (Andersson et al., 2007; Fantozzi et al., 2013; Gårdfeldt et al., 2003), and the atmospheric deposition of HgT is set to 67.5 Mg yr−1, 2 % of which is MMHg (Cossa et al., 2022, 2017). The Hg concentrations in the Black Sea outflow at the Dardanelles Strait are set at 0.2 pM in agreement with the observations from the GEOTRACES GA04 cruise and modeled fluxes (Rosati et al., 2018).

2.3 Model sensitivity to key parameters and processes

A set of sensitivity simulations was implemented to investigate model uncertainty and improve the fit against experimental observations, initially showing substantial underestimation of MeHg concentration maxima. The sensitivity simulations aimed at assessing the impacts of variations of the values of the sinking velocity of organic detritus and associated Hg species (Sect. 2.3.1) and of the coefficient for Hg methylation rates (Sect. 2.3.2), as well as the loading of Hg from rivers (Sect. 2.3.3), which is a process poorly characterized by field observations (Liu et al., 2021). All the analyses of model results focus on the same simulation, which adopts a sinking velocity of 10 m d−1, a methylation rate constant higher than the one used in global modeling (Zhang et al., 2020), and accounts only for dissolved input of Hg from rivers.

2.3.1 Sensitivity simulations for POC sinking velocity

A sensitivity analysis of POC sinking velocity was carried out, hypothesizing that a change in this parameter would impact the vertical distribution of all Hg species due to a change in the distribution and remineralization of POC and in the transport of particulate Hg and MeHg. We also speculate that a deeper occurrence of POC remineralization, leading to deeper MeHg maximum, could result in a decreased amount of MeHg photochemically degraded and thus in higher MeHg concentrations. The current version of the BFM model includes only one class of non-living organic detritus (POC) that is a sink for all the plankton state variables (Vichi et al., 2015), resulting in a homogenous POC sinking velocity representative of an ecosystem-averaged dynamic. In the field, organic particles and aggregates that originate from a continuous size spectrum of planktonic cells of various composition show a wide range of sinking speeds (Cael et al., 2021). The reference model simulation was run by adopting the default POC sinking velocity of 3 m d−1 previously used to model the biogeochemistry of the Mediterranean Sea (Lazzari et al., 2021, 2012). Here, two additional sensitivity simulations were run by setting the POC sinking velocity at 10 and 20 m d−1, assuming community-averaged cell diameters of 35 and 50 µm, respectively.

2.3.2 Sensitivity simulations for Hg methylation

Hg methylation was parameterized in the model by scaling the flux of POC remineralization for a constant xmet (Table S2), adopting the same value used in a global ocean model (Zhang et al., 2020), due to the unavailability of more specific data. We explored the sensitivity of modeled MeHg concentrations to a 3-fold increase in xmet, as well as to an alternative parameterization that included both POC and DOC remineralization flux (Eq. S10).

2.3.3 Sensitivity simulation for river Hg load

The largest Mediterranean rivers are characterized by proximal-accumulation-dominated dispersal system, i.e., deltaic systems with fast and substantial (> 50 %–90 %) sediment accumulation in the proximity of their mouth (Walsh and Nittrouer, 2009). Consistently with this picture, a recent budget for the Mediterranean Sea estimated that the magnitude of total Hg input from rivers (6 Mg yr−1) is comparable to the Hg flux to the shelf sediments (6.8 Mg yr−1), suggesting that most of the Hg associated with riverine particles settles before reaching the open ocean (Cossa et al., 2022). We, therefore, included only the dissolved load of Hg and MMHg in the reference simulation and explored the uncertainty related to this choice with a sensitivity simulation accounting for the total load Hg species, i.e., including both the particulate and dissolved inputs of HgII and MMHg. Limited data are available to characterize the spatial and temporal variability of these sources. The most robust observations are from a time series (2008–2010) for the Rhône River (Cossa et al., 2017) indicating 0.85 ± 0.45 nmol g−1 of HgT and 0.017 ± 0.008 nmol g−1 of MMHg. Assessments in the Ebro River, impacted by an industrially polluted water reservoir, and in the Po River, which drain the largest Italian industrial and agricultural areas, found respectively 4.9 and 1.15 nmol g−1 of particulate HgT at the river mouth, also highlighting the role of episodic stormy events in transporting large amounts of particles and associated pollutants (Palanques et al., 2020; Vignati et al., 2003). Further, at the mouth of the Soča/Isonzo River, impacted by upstream Hg mining, concentrations of unfiltered Hg higher than 200 pM have been reported (Hines et al., 2000). Hg loads for the Mediterranean tributary were estimated based on this information as synthesized in Table S6, also qualitatively referring to the spatial distribution of Hg stock and erosion fluxes in European topsoils (Panagos et al., 2021) and other information on anthropogenic activities (Ahmed et al., 2018). In the reference simulation, including only the dissolved input, the loading estimate is 0.28 Mg yr−1 for total dissolved Hg (and 0.008 Mg yr−1 for dissolved MMHg), while in the simulation with total Hg inputs the loading is about 4.6 Mg yr−1, including 0.27 Mg yr−1 of MMHg. Despite the relatively high mean concentrations (Table S6) assumed, the estimated HgT load is slightly lower than the estimate (6 Mg yr−1) from Cossa et al. (2022), most likely because the water discharge used in that work (> 350 km3 yr−1) is higher than the total discharge used in this model setup (> 291 km3 yr−1). On the other hand, the load to the Mediterranean Sea from solely European rivers (2.9 Mg yr−1) estimated from Hg concentrations in topsoils and water erosion rates from river catchments (Panagos et al., 2021) is slightly lower than the load from European rivers in our model setup (3.25 Mg yr−1). Finally, all these estimates are lower than the 18 Mg yr−1 computed by Liu et al. (2021) for the Mediterranean Sea through a global river model explicitly considering water discharge, suspended solids, and Hg loadings. Based on these calculations, Liu et al. (2021) proposed a global river load (1000 Mg yr−1) much higher than the estimate from the 2018 global mercury assessment (300 Mg yr−1) (UNEP, 2019). They also estimated that 70 % of the load is buried in coastal and shelf sediment.

3 Results and discussion

3.1 Model sensitivity and uncertainty

The sensitivity analysis showed that the sinking rate of organic detritus (POC) is an important control on the vertical distribution of MeHg (i.e., the shape of the vertical profile) along the water column but has little effect on the values of MeHg concentration maxima and the distribution of inorganic Hg species (Supplement Sect. S1.1 and Figs. S1 and S2). On the other hand, the increase in sinking velocity from 3 to 10 m d−1 improved the agreement of modeled seasonal fluxes of POC at 200 m depth with available observations (Ramondenc et al., 2016) (Supplement Sect. S1.1 and Fig. S3). The sensitivity analysis also showed that an increase in the coefficient for Hg methylation (xmet) results in a higher agreement with observed MeHg concentrations (Cossa et al., 2022, 2009), while the inclusion of DOC remineralization flux results in an overestimation of surficial MeHg concentrations (Supplement Sect. S1.2 and Fig. S3).

The inclusion of higher riverine HgII, MMHg, and POC inputs does not improve the agreement between modeled and observed concentrations of MeHg (Supplement Sect. S1.3 and Figs. S5–S10), suggesting that the underestimation of MeHg maxima in subsurface waters is not due to an underestimation of watershed sources. These results also corroborate the idea that a substantial fraction of particulate Hg and MeHg carried by rivers settles in coastal areas along with POC (Fig. S5), as estimated in a recent budget based on field observations (Cossa et al, 2022). Nonetheless, large uncertainties in the spatial and temporal variability of riverine loadings and coastal processes remain to be addressed both at the basin and global scale (Cossarini et al., 2021; Liu et al., 2021).

The underestimation of MeHg in the model does not appear to be related to inorganic Hg availability, since the concentrations of HgT are within the experimental uncertainty of the data (Fig. S11).

Modeled fluxes of HgT and Hg0 compare well with other budget estimates (Cossa et al., 2022; Žagar et al., 2007). The net exchange of gaseous Hg0 with the atmosphere indicates 44.8 Mg yr−1, which is close to a previous estimate (50 Mg yr−1) based on aquatic biogeochemical modeling (Žagar et al., 2007) and slightly lower than the flux adopted in the budget from Cossa and coauthors (67.5 Mg yr−1). Considering both Hg0 exchange and atmospheric deposition of HgII (Sect. 2.2.6 and Fig. 3), the net exchange of Hg species ( 7 Mg yr−1) is here significantly lower than the estimate (29.8 Mg yr−1) from Cossa et al. (2022). However, while the latter estimate is based on an atmospheric model assuming constant ocean surface Hg0 concentrations (Gencarelli et al., 2014), here, the flux is estimated accounting for 3D dynamics in the marine system, including sequestration of Hg species from surface waters due to the biological carbon pump. Moreover, recent evidence from isotopes studies pointed out that the ocean might be a net sink for atmospheric Hg rather than a net source (Jiskra et al., 2021).

The modeled export of HgT from the Alboran Sea to the Atlantic Ocean is 9.3 Mg yr−1, close to the value (9.7 Mg yr−1) estimated by Cossa et al. (2022). On the other hand, our estimate (2.5 Mg yr−1) of MeHg outflow is significantly lower than the 5.6 Mg yr−1 of the budget by Cossa et al. (2022), even assuming an increased value for xmet and, consequently, of subsurface MeHg.

All in all, a calibration of the parameter xmet appears appropriate based on the current understanding, also considering that the initial guess for this parameter comes from a global ocean application (Zhang et al., 2020). Indeed, such a global model study underestimated the observed MeHg concentrations in the Mediterranean Sea, while reproducing with good agreement, or overestimating, observations from various cruises in the North Atlantic and Pacific Ocean, pointing out the need to better resolve spatial and temporal variability of biological processes and linked Hg dynamics (Zhang et al., 2020). A comprehensive validation of the OGSTM–BFM biogeochemical variables (Cossarini et al., 2021) highlighted good model accuracy in surface water but higher mismatches with observed oxygen and nutrient concentrations in the mesopelagic zone, attributed to uncertainties in remineralization and deep ventilation processes, which also affect methylation dynamics. In spite of uncertainties remaining to be addressed, the OGSTM–BFM–Hg model is able to reproduce the spatial gradients observed in the Mediterranean Sea, showing the improvements in capability to simulate Hg dynamics (Žagar et al., 2007, 2014) through the coupling with key biological processes, e.g., the carbon pump and microbial loop (Bowman et al., 2020; Cossa et al., 2017, 2009; Heimbürger et al., 2010; Monperrus et al., 2007; Munson et al., 2018; Sunderland et al., 2009).

3.2 MeHg distribution in the water of the Mediterranean Sea

The zonal and vertical gradients of MeHg concentrations reproduced by the OGSTM–BFM–Hg model show the highest concentrations in subsurface waters of the most productive subbasins of the Mediterranean Sea (Fig. 4), consistently with the available observations (Cossa et al., 2009). Modeled vertical profiles of MeHg concentrations (Fig. 5) are in good agreement with the observations in the Southern Adriatic Sea (Sad, Fig. 5a) and are in the lower range of the observations for the Northwestern Mediterranean (Nwm, Fig. 5b) and Ionian (Ion, Fig. 5c) subbasins. Seasonality affects the modeled distribution of MeHg (Figs. 4 and 5) with more marked effects in the Sad subbasin (Fig. 5a). The strong winter convection occurring in the Sad in 2017 caused remixing in the water column (Mihanović et al., 2021), leading to the disruption of the subsurface MeHg maximum from January to March and an increase in MeHg surface water concentrations, in spite of net demethylation prevailing at all depths during winter months (Fig. 6). A slighter increase in surface water MeHg concentrations during winter months is visible also for the other subbasins (Figs. 4, 5b, and c), but, in the absence of strong winter convection the variation is small, and the model predicts that the subsurface MeHg maximum is a permanent feature throughout the year. The Nwm is known to be an area of winter convection; however, other authors reported a recent decline of deep convection for this area that was weak in 2014 and never occurred in the period 2015–2017 (Margirier et al., 2020). Based on these model results, we speculate that in the long run the reduction of winter mixing can cause MeHg accumulation in the biologically active zone enhancing the exposure of marine organisms.

Figure 4Spatial distribution of monthly averaged MeHg concentrations (pM) in the Mediterranean Sea at different depths in the water column (25, 90, 175, 270, and 710 m depth) simulated with the OGSTM–BFM–Hg model.


In the Sad subbasin, a progressive buildup of MeHg subsurface maxima (up to 0.18 pM) is predicted from April to September, driven by an increase in primary production triggering higher POC remineralization and in turn higher Hg methylation rate constants kmet (Fig. S4). However, while kmet values are maxima in April and net methylation fluxes (Fig. 6) are maxima in June (0.45 pmol m−3 d−1), MeHg concentrations are maxima in September and remain stable until October (Fig. 5), highlighting the importance of water stratification for the formation of the subsurface maxima. Figure 6 also shows that net demethylation prevails throughout the year in deep waters (except for weak positive methylation in Sad and Tyr during summer) and for most of the year in surface waters, suggesting that the presence of MeHg in surface and deep water depends on its formation in the intermediate waters and diapycnal mixing.

Modeled vertical profiles of MeHg in the oligotrophic Ion subbasin (Fig. 5c) have lower maxima (up to 0.16 pM) and a smoother concentration gradient compared to the more productive waters of the Nwm (up to 0.24 pM). The west–east gradient is sustained by the highest methylation rates (Fig. S12) in the western Mediterranean, driven by higher primary productivity and remineralization of organic detritus, and it is further reinforced by high rates of both photochemical and biological MeHg degradation in the eastern subbasins (Fig. S12), due to the highest water temperature, irradiance, and deeper penetration of shortwave radiation caused by the low plankton biomass.

Figure 5Vertical profiles of monthly averaged MeHg (MMHg + DMHg) concentrations in the water column reproduced by the OGSTM–BFM–Hg model for 2017 (colored solid lines), compared with observations for 2004–2005 (dotted lines with crosses) available from the literature (Cossa et al., 2009) in three subbasins (Fig. 1): (a) Southern Adriatic Sea (Sad), (b) Northwestern Mediterranean (Nwm), and (c) Ionian Sea (Ion2). The right panels show the entire water column, and the left panels show a detail of modeled concentrations in the top 600 m of the same subbasins.


Figure 6Monthly evolution of the net methylation fluxes (pmol m−3 d−1) for each subbasin (colored markers) of the Mediterranean Sea (black lines) depth-integrated for (a) surface water (0–100 m depth), (b) intermediate water (100–600 m depth), and (c) deep water (> 600 m depth). Net methylation is calculated as the difference between modeled Hg methylation fluxes and demethylation fluxes (both dark and photochemical).


3.3 Modeled MMHg in phytoplankton and zooplankton

The spatial and temporal distributions of modeled MMHg in phytoplankton (MMHgphy) and zooplankton (MMHgzoo) of the Mediterranean Sea (Fig. 7) are overall comparable to the concentrations of particulate MMHg (0.7 ± 1.3 pmol m−3) measured in the Atlantic Ocean (Bowman et al., 2015). The model reproduced a seasonal increase in MMHgphy (Fig. 7a) during late spring and summer that is driven by blooms of picophytoplankton (Fig. S13) combined with an increase in MMHg availability in surface water (Fig. S14). The amount of MMHg bioaccumulated by picophytoplankton (MMHgphy,P3, Fig. S15) is high (up to 0.6 pmol(Hg) m(w)-3) compared to the bioaccumulation by other phytoplankton PFTs (up to 0.02 pmol(Hg) m(w)-3, Figs. S16–S18) due to the model assumption of increasing bioconcentration capacity with decreasing phytoplankton cell size (Sect. 2.3). However, not all the picophytoplankton blooms result in high values of MMHgphy: the strongest blooms of picophytoplankton, occurring in March and April at 25 m depth (Fig. S13), have a marginal impact on bioaccumulation, especially in the easternmost subbasins, due to the very low concentrations of MMHg in the most surficial waters (Fig. S14). On the other hand, much weaker blooms of picophytoplankton occurring from May to August slightly deeper in the water column (35 and 45 m depth) lead to a maximum in phytoplankton MMHg concentrations (Fig. 7a) because of higher water MMHg concentrations at these depths during these months (Fig. S14).

Figure 7Spatial distribution of monthly averaged (a) MMHgphy and (b) MMHgzoo (pmol(Hg) m(w)-3) in the Mediterranean Sea at different water depths (25, 35, and 45 m) simulated with the OGSTM–BFM–Hg model for 2017. The plots show the sum of all variables of MMHg in phytoplankton (Figs. S7–S10) and zooplankton (Figs. S11–S14).


The average MMHgphy,P3 is lower in the Alb subbasin than in the Sad subbasin (Table 1), the latter being characterized by higher availability of MMHg at 35–45 m depth and higher biomasses of picophytoplankton during summer (Figs. S13 and S14). The hypothesis of the cell size effect on bioconcentration is the best synthesis of available knowledge (Schartup et al., 2018); however, the experimental study from which this relationship was derived also showed deviations from this pattern and temperature-induced effects, possibly due to phagocytosis of Hg-OM compounds by mixotroph organisms such as dinoflagellates (Lee and Fisher, 2016). One of the few field studies assessing the size effect on Hg bioaccumulation in plankton (Gosnell and Mason, 2015) found an increase in MeHg % but not in MeHg content for smaller plankton size (< 5 µm); however, the size classes analyzed by that study do not overlap well with the functional groups of the OGSTM–BFM–Hg model. More field, in vitro, and modeling studies are needed to disentangle the dynamics underlying key processes for Hg bioaccumulation.

Table 1Indicators for MMHg bioaccumulation compared for the Alb and Sad subbasins (Fig. 3). The average (± SD) values of the plankton MMHg content (MMHgPFT, pmol m−3), the MMHg quota (QMMHg, PFT, ng gw.w.-1), and the bioconcentration factors (log(BCFPFT)) are shown for each phytoplankton and zooplankton PFT.

Download Print Version | Download XLSX

The highest concentrations of MMHgzoo (Fig. 7b) are predicted from June to September at 35 m depth, showing a delay of about 1 month with respect to the temporal evolution of MMHgphy due to the combined impacts of plankton phenology and trophic interactions. Heterotrophic nanoflagellates bioaccumulate much more MMHg (up to 0.5 pmol m−3, Fig. S19) than other zooplankton PFTs (up to 0.1 pmol m−3, Figs. S20–S22) because they feed on picophytoplankton (Fig. 8), which is more enriched in MeHg than other PFTs and is an important driver for the distribution of MMHgzoo. However, in spite of heterotrophic nanoflagellates being always abundant from March to October (Fig. S23) and most abundant in the eastern subbasins, the highest MMHg bioaccumulation occurs when the increase in the nanoflagellates population is coincident with an enrichment of MMHg in phytoplankton (i.e., in western subbasins, as well as in the Adriatic and part of the Ionian Sea during summer).

Figure 8Conceptual model of MMHg bioaccumulation in the plankton food web of the OGSTM–BFM–Hg model, subdivided into a lower “herbivorous” (red arrows) and upper “omnivorous” (blue arrows) food webs. The squares with bold numbers show estimated values of the bioconcentration factors (log(BCF)) for each PFT, calculated for the Sad subbasin from spatially averaged model output at monthly resolution. The green circles indicate the four phytoplankton PFTs (P1, P2, P3, and P4, representative respectively of diatoms, autotrophic nanoflagellates, picophytoplankton, and large phytoplankton), the light blue circles indicate the four zooplankton PFTs (Z3, Z4, Z5, and Z6, representative respectively of carnivorous mesozooplankton, omnivorous mesozooplankton, microzooplankton, and heterotrophic nanoflagellates), the purple circle indicates bacteria (B1, representative of all heterotrophic bacteria), and the grey circles indicate particulate organic carbon (POC, representative of all organic detritus) and semi-refractory dissolved organic matter (DOCsr). The solid arrows show the preferential trophic interactions, the dashed arrows are for secondary trophic interactions (diet switch if the favorite preys are not available), and the circle dashed arrows are for cannibalization within the same zooplankton group. The dash and dotted arrows indicate passive uptake of MMHg by phytoplankton. The dotted arrows indicate the processes of transformation related to carbon recycling by bacteria (remPOC) and of Hg methylation, whose kinetics are coupled through the parameter kmet, converting HgII into MMHg.


The mean MMHg quota (QMMHg,PFT) in phytoplankton of the Mediterranean Sea is about 3.0 ng gw.w.-1 for picophytoplankton, 0.3 ng gw.w.-1 for autotrophic nanoflagellates, 0.03 ng gw.w.-1 for diatoms, and 0.04 ng gw.w.-1for large phytoplankton. The spatial variability is highlighted by comparing the quotas estimated for the Alb and Sad subbasins (Table 1) that are in the lower and upper range of the values for the Mediterranean Sea. Modeled values are comparable with the few observations available for open waters of the central Pacific Ocean (0.1–4 ng gw.w.-1) (Gosnell and Mason, 2015) and of the northwest Atlantic Ocean (0.15 ± 0.06 ng gw.w.-1) (Hammerschmidt et al., 2013), as well as with data from various coastal and shelf ecosystems (Harding et al., 2018). The speciation of Hg in phytoplankton groups of the Mediterranean Sea has never been assessed (Cinnirella et al., 2019). The MMHg quota was about 0.14 ± 0.1 ng gw.w.-1 in seston of diameter 80–200 µm sampled in the Gulf of Lion during spring and autumn 2004–2006 (Cossa et al., 2012).

The mean QMMHg,PFT values in zooplankton are 0.8 ng gw.w.-1 for heterotrophic nanoflagellates, 0.3 ng gw.w.-1 for microzooplankton, 0.6 ng gw.w.-1 for omnivorous mesozooplankton, and 0.9 ng gw.w.-1 for carnivorous mesozooplankton. The MMHg quota in carnivorous zooplankton (QMMHg,Z3) estimated for the Sad subbasin is much higher (3.60 ± 2.2) than the average value due to the high MMHg accumulation in the whole food web (Table 1). However, since biomasses of carnivorous zooplankton in the Sad are very low (Fig. S16), the overall content of MMHg in this group (MMHgzoo,Z3) is lower (1.0 × 10−3± 2.7 × 10−4) than in the Alb subbasin (1.9 × 10−2± 1.6 × 10−3).

The modeled zooplankton quotas compare well with the observations of MMHg in mesozooplankton of the Gulf of Lion, which was 0.52 ± 0.39 ng gw.w.-1 for the size range 200–500 µm, 0.57 ± 0.39 ng gw.w.-1 for the size range 500–1000 µm, and 1.8 ± 1.3 ng gw.w.-1 for zooplankton > 200 µm (Cossa, 2012). MMHg measured in the 1990s in zooplankton in the Gulf of Trieste (Nad, Fig. 1) Hg hotspot was 2.1 ± 2.3 ng gw.w.-1 (Horvat et al., 1999, 2014), consistently with the fact that the Gulf of Trieste is an upper bound for the Mediterranean, being subjected to high loadings of legacy Hg from the Soča/Isonzo river (Hines et al., 2000) and coastal lagoons (Melaku Canu et al., 2015; Rosati et al., 2020).

3.4 Bioconcentration factors and trophic magnification factors

Modeled log(BCFs) values for phytoplankton PFTs of the Mediterranean Sea range 3.7–6.5, which is slightly lower than the range (4.3–6.8) reported for phytoplankton from the central Pacific (Gosnell and Mason, 2015) and is higher than the ranges for Long Island Sound (2.6–5.5) (Gosnell et al., 2017) and other coastal and shelf sites of the world (3.1–4.4) (Harding et al., 2018). Higher BCFs in open waters than in coastal areas are consistent with the enhanced Hg bioavailability in the open sea inferred by other authors (Gosnell and Mason, 2015; Schartup et al., 2015), which in the model is approximated through the inverse relationship between MMHg uptake and DOC (Eq. 10). Modeled zooplankton log(BCFs) values in the Mediterranean Sea span from 4.0–6.4, comparable to the range for mesozooplankton (4.0–6.5) in the Pacific Ocean (Gosnell and Mason, 2015).

The trophic interactions among PFTs simulated by the model can be subdivided into a lower herbivorous part of the food web and an upper omnivorous part of the food web (Fig. 8).

In the lower part of the food web, heterotrophic nanoflagellates (2–20 µm) graze on bacteria (< 1 µm) and picophytoplankton (0.2–2 µm), while microzooplankton (20–200 µm) grazes on heterotrophic and autotrophic nanoflagellates (2–20 µm), diatoms (200–20 µm), and occasionally on other PFTs. In the upper part of the food web, the omnivorous mesozooplankton (200–2000 µm) preys on microzooplankton and all the phytoplankton PFTs > 2 µm, while carnivorous mesozooplankton (200–2000 µm) preys on omnivorous mesozooplankton and large phytoplankton (> 100 µm). The highest BCFs (and QMMHg,PFT) are predicted for picophytoplankton and carnivorous zooplankton (Table 1 and Fig. 8). Heterotrophic nanoflagellates display lower bioaccumulation than picophytoplankton because they also feed on bacteria and autotrophic nanoflagellates that dilute the MMHg intake (Fig. 8). Microzooplankton and omnivorous zooplankton, which have an intermediate position in the food web, have bioaccumulation levels comparable to those of heterotrophic nanoflagellates. The effect of biomagnification is visible only in carnivorous zooplankton, having much higher BCFs than their prey (Table 1). These dynamics are consistent with a global model taking into account a more complex plankton food web (366 PFTs) that found higher MMHg levels in phytoplankton than in their predators and biomagnification in the highest trophic levels with high food intake (Wu et al., 2021).

Figure 9Trophic magnification factor (TMF, unitless) for the (a) omnivorous and (b) carnivorous plankton food webs (Fig. 8) for each subbasin of the Mediterranean Sea (Fig. 1), calculated by integrating over 0–100 m depth the monthly outputs of the OGSTM–BFM–Hg model.


Trophic magnification factors (TMFs) were calculated for each subbasin of the Mediterranean Sea for the lower and upper parts of the food web. Low TMFs (range 0.05–0.38) were estimated for the lower food web, indicating the absence of biomagnification (Fig. 9a), and higher TMFs (range 0.15–2.60) were estimated for the upper food web (Fig. 9b), in agreement with previous modeling studies suggesting that MMHg biomagnification between small zooplankton groups and phytoplankton is unlikely to happen (Wu et al., 2021, 2020; Zhang et al., 2020). Biomagnification (TMF > 1) in the upper food web (Fig. 9b) occurs only in subbasins with relatively abundant biomasses of carnivorous zooplankton (Fig. S24), namely in the Alb and Swm subbasins. The TMF for the Alb subbasin (2.4 ± 0.28) is lower than but comparable to the biomagnification from microseston to zooplankton (TMF = 4) estimated from field data for the northwest Atlantic Ocean (Hammerschmidt et al., 2013). The other Mediterranean subbasins, due to the absence or very low abundance of carnivorous mesozooplankton (< 0.05 mg(C) m−3), have a shorter food web and lower TMFs (< 1). Although oligotrophy tends to favor longer food webs enhancing biomagnification, in extremely oligotrophic ecosystems the very low primary production can limit the presence of carnivorous zooplankton populations, lowering the bioaccumulation potential (Wu et al., 2021). The TMFs < 1 for the upper part of the food web are comparable to the estimates made with a global model for the most oligotrophic ocean regions such as the South Pacific Gyre, the North Atlantic Gyre, and the South Atlantic Gyre (Wu et al., 2021).

4 Conclusions

We developed and released the OGSTM–BFM–Hg model, a numerical state-of-the-art model tracking the biogeochemistry of the main marine Hg species (HgII, Hg0, MMHg, and DMHg) fully coupled to the biogeochemistry of plankton and nutrients. The model, applied to a high-resolution 3D domain of the Mediterranean Sea, successfully captured the zonal and vertical variability of MeHg concentrations, reproducing a decreasing MeHg gradient from west to east that is supported by the available observations. Sensitivity analysis showed that the agreement with observed MeHg concentrations and POC fluxes was improved by adopting an intermediate sinking speed for organic detritus. Modeled MeHg concentrations also improved by increasing the coefficient constant for the Hg methylation rate, while the inclusion of the dissolved organic carbon remineralization in the equation for Hg methylation increases the mismatch between the model and observations.

Sensitivity simulations on river Hg loadings corroborate the idea that in the Mediterranean Sea most of the inputs of particulate Hg (and POC) from rivers sink in coastal and shelf areas and have a limited impact on open-sea dynamics. However, as pointed out in a recent reassessment of global river inputs, there is a need to better assess the spatial and temporal variability of riverine loadings, especially in relation to extreme events and their increasing occurrence driven by climate changes. The net Hg exchange from Mediterranean waters to the atmosphere is lower than in a previous budget but in line with findings from a recent Hg isotopes study.

Model results highlight that summer stratification of the water column is an important process for the buildup of a subsurficial maximum of MeHg concentrations, while the occurrence of deep convection in winter results in substantial redistribution of MeHg smoothing the vertical profiles. A decrease in winter convection events linked to increasing water temperature, which has already been observed in recent years in the Mediterranean Sea, seems to limit MeHg redistribution in the water column, causing higher water MeHg concentrations in the biologically active zone and, in turn, higher plankton exposure. In fact, spatial and temporal variability of plankton bioaccumulation in the model is controlled by plankton phenology and by the availability of MMHg at the depths at which plankton blooms occur. The biomagnification potential is stronger in productive areas of the Mediterranean Sea characterized by high biodiversity and longer food web length with a significant presence of carnivorous zooplankton, such as the Alb and Swm subbasins.

Code availability

The OGSTM–BFM–Hg model code is publicly available in the Zenodo repository at (Rosati et al., 2022a).

Data availability

The data used for model implementation and validation are available upon request.


The supplement related to this article is available online at:

Author contributions

The conceptualisation and investigation were carried out by GR, DC, PL, and CS. The software and methodology was developed by GR, DC, PL, and CS. The formal analysis, validation, and visualization were performed by GR and PL. Funding acquisition and project administration were carried out by DC and CS. Supervision and resources were provided by DC and CS. Writing of the original draft was performed by GR. Review and editing were performed by GR, DC, PL, and CS.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


We thank the reviewers for their constructive comments that improved the quality of the work.

Financial support

The research has been partially supported by the PRIN project ICCC (Impacts of climate change on the biogeochemistry of contaminants in the Mediterranean Sea), funded by the Ministero dell'Istruzione, dell'Università e della Ricerca (grant no. 2017ZRPNKX_001), and by the Interreg MED Strategic Project SHAREMED, co-financed by the European Regional Development Fund under the Funding Programme Interreg MED 2014–2020.

Review statement

This paper was edited by Gwenaël Abril and reviewed by Yanxu Zhang and Gwenaël Abril.


Ahmed, A. M. A., Purwanto, P., and Sunoko, H. R.: Consequences of Mercury Used by Artisanal and Small-Scale Gold Mining Processes a Case of River Nile State Sudan, J. Ecol. Eng., 20, 106–115, 2019. 

Alava, J. J., Cisneros-Montemayor, A. M., Sumaila, U. R., and Cheung, W. W. L.: Projected amplification of food web bioaccumulation of MeHg and PCBs under climate change in the Northeastern Pacific, Sci. Rep., 8, 1–12,, 2018. 

Amos, H. M., Sonke, J. E., Obrist, D., Robins, N., Hagan, N., Horowitz, H. M., Mason, R. P., Witt, M., Hedgecock, I. M., Corbitt, E. S., and Sunderland, E. M.: Observational and modelling constraints on global anthropogenic enrichment of mercury, Environ. Sci. Technol., 49, 4036–4047,, 2015. 

An, J., Zhang, L., Lu, X., Pelletier, D. A., Pierce, E. M., Johs, A., Parks, J. M., and Gu, B.: Mercury Uptake by Desulfovibrio desulfuricans ND132: Passive or Active?, Environ. Sci. Technol., 53, 6264–6272,, 2019. 

Andersson, M. E., Gårdfeldt, K., Wängberg, I., Sprovieri, F., Pirrone, N., and Lindqvist, O.: Reprint of “Seasonal and daily variation of mercury evasion at coastal and off shore sites from the Mediterranean Sea”, Mar. Chem., 107, 104–116,, 2007. 

Baya, P. A., Gosselin, M., Lehnherr, I., St. Louis, V. L., and Hintelmann, H.: Determination of monomethylmercury and dimethylmercury in the arctic marine boundary layer, Environ. Sci. Technol., 49, 223–232,, 2015. 

Blum, J. D., Popp, B. N., Drazen, J. C., Anela Choy, C., and Johnson, M. W.: Methylmercury production below the mixed layer in the North Pacific Ocean, Nat. Geosci., 6, 879–884,, 2013. 

Bowman, K. L., Hammerschmidt, C. R., Lamborg, C. H., and Swarr, G.: Mercury in the North Atlantic Ocean: The U.S. GEOTRACES zonal and meridional sections, Deep-Sea Res. Pt. II, 116, 251–261,, 2015. 

Bowman, K. L., Hammerschmidt, C. R., Lamborg, C. H., Swarr, G. J., and Agather, A. M.: Distribution of mercury species across a zonal section of the eastern tropical South Pacific Ocean (U.S. GEOTRACES GP16), Mar. Chem., 186, 156–166,, 2016. 

Bowman, K. L., Lamborg, C. H., and Agather, A. M.: A global perspective on mercury cycling in the ocean, Sci. Total Environ., 710, 136166,, 2020. 

Bosc, E., Bricaud, A., and Antoine, D.: Seasonal and interannual variability in algal biomass and primary production in the Mediterranean Sea, as derived from 4 years of SeaWiFS observations, Global Biogeochem. Cy., 18, 1–17,, 2004. 

Bricaud, A., Bosc, E., and Antoine, D.: Algal biomass and sea surface temperature in the Mediterranean Basin Intercomparison of data from various satellite sensors, and implications for primary production estimates, Remote Sens. Environ., 81, 163–178,, 2002. 

Cael, B. B., Cavan, E. L., and Britten, G. L.: Reconciling the Size-Dependence of Marine Particle Sinking Speed, Geophys. Res. Lett., 48, 1–11,, 2021. 

Canu, D. and Rosati, G.: Long-term scenarios of mercury budgeting and exports for a Mediterranean hot spot (Marano-Grado Lagoon, Adriatic Sea), Estuar. Coast. Shelf Sci., 198, 518–528,, 2017. 

Canu, D., Ghermandi, A., Nunes, P. A. L. D., Lazzari, P., Cossarini, G., and Solidoro, C.: Estimating the value of carbon sequestration ecosystem services in the mediterranean sea: An ecological economics approach, Glob. Environ. Change, 32, 87–95,, 2015. 

Canu, D. M., Rosati, G., and Solidoro, C.: Mercury Budget and Scenario Analysis for the Marano-Grado Lagoon, Using Modelling and Observations, Proc., 30, 19,, 2019. 

Choe, K.-Y. and Gill, G. A.: Distribution of particulate, colloidal, and dissolved mercury in San Francisco Bay estuary, 2. Monomethyl mercury, Limnol. Oceanogr., 48, 1547–1556,, 2003. 

Choe, K.-Y., Gill, G. A., and Lehman, R.: Distribution of particulate, colloidal, and dissolved mercury in San Francisco Bay estuary, 1. Total mercury, Limnol. Oceanogr., 48, 1535–1546,, 2003. 

Cinnirella, S., Bruno, D. E., Pirrone, N., Horvat, M., Živković, I., Evers, D. C., Johnson, S., and Sunderland, E. M.: Mercury concentrations in biota in the Mediterranean Sea, a compilation of 40 years of surveys, Sci. Data, 6, 205,, 2019. 

Coale, K. H., Heim, W. A., Negrey, J., Weiss-Penzias, P., Fernandez, D., Olson, A., Chiswell, H., Byington, A., Bonnema, A., Martenuk, S., Newman, A., Beebe, C., and Till, C.: The distribution and speciation of mercury in the California current: Implications for mercury transport via fog to land, Deep-Sea Res. Pt. II, 151, 77–88,, 2018. 

Conaway, C. H., Black, F. J., Gault-Ringold, M., Pennington, J. T., Chavez, F. P., and Flegal, A. R.: Dimethylmercury in coastal upwelling waters, Monterey Bay, California, Environ. Sci. Technol., 43, 1305–1309,, 2009. 

Cossa, D. and Coquery, M.: The Mediterranean Mercury Anomaly, a Geochemical or a Biologocal Issue, in: The Mediterranean Sea, edited by: Saliot, A., Springer Berlin Heidelberg, Berlin, Heidelberg, 177–208, ISBN 978-3-540-31492-9, 2005. 

Cossa, D., Averty, B., and Pirrone, N.: The origin of methylmercury in open Mediterranean waters, Limnol. Oceanogr., 54, 837–844,, 2009. 

Cossa, D., Harmelin-Vivien, M., Mellon-Duval, C., Loizeau, V., Averty, B., Crochet, S., Chou, L., and Cadiou, J.-F.: Influences of bioavailability, trophic position, and growth on methylmercury in hakes (Merluccius merluccius) from Northwestern Mediterranean and Northeastern Atlantic, Environ. Sci. Technol., 46, 4885–93,, 2012. 

Cossa, D., Durrieu de Madron, X., Schäfer, J., Lanceleur, L., Guédron, S., Buscail, R., Thomas, B., Castelle, S., and Naudin, J.-J.: The open sea as the main source of methylmercury in the water column of the Gulf of Lions (Northwestern Mediterranean margin), Geochim. Cosmochim. Ac., 199, 222–237,, 2017. 

Cossa, D., Knoery, J., Bănaru, D., Harmelin-Vivien, M., Sonke, J. E., Hedgecock, I. M., Bravo, A. G., Rosati, G., Canu, D., Horvat, M., Sprovieri, F., Pirrone, N., and Heimbürger-Boavida, L.-E.: Mediterranean Mercury Assessment 2022: An Updated Budget, Health Consequences, and Research Perspectives, Environ. Sci. Technol., 56, 3840–3862,, 2022. 

Cossarini, G., Lazzari, P., and Solidoro, C.: Spatiotemporal variability of alkalinity in the Mediterranean Sea, Biogeosciences, 12, 1647–1658,, 2015. 

Cossarini, G., Feudale, L., Teruzzi, A., Bolzon, G., Coidessa, G., Solidoro, C., Di Biagio, V., Amadio, C., Lazzari, P., Brosich, A., and Salon, S.: High-Resolution Reanalysis of the Mediterranean Sea Biogeochemistry (1999–2019), Front. Mar. Sci., 8, 1–21,, 2021. 

Crise, A., Allen, J. I., Baretta, J., Crispi, G., Mosetti, R., and Solidoro, C.: The Mediterranean pelagic ecosystem response to physical forcing, Prog. Oceanogr., 44, 219–243,, 1999. 

Crispi, G., Mosetti, R., Solidoro, C., and Crise, A.: Nutrients cycling in Mediterranean basins: The role of the biological pump in the trophic regime, Ecol. Modell., 138, 101–114,, 2001. 

De Simone, F., Gencarelli, C. N., Hedgecock, I. M., and Pirrone, N.: Global atmospheric cycle of mercury: a model study on the impact of oxidation mechanisms, Environ. Sci. Pollut. Res., 21, 4110–4123,, 2014. 

Di Biagio, V., Cossarini, G., Salon, S., Lazzari, P., Querin, S., Sannino, G., and Solidoro, C.: Temporal scales of variability in the Mediterranean Sea ecosystem: Insight from a coupled model, J. Mar. Syst., 197, 103176,, 2019. 

Di Biagio, V., Cossarini, G., Salon, S., and Solidoro, C.: Extreme event waves in marine ecosystems: An application to Mediterranean Sea surface chlorophyll, Biogeosciences, 17, 5967–5988,, 2020. 

D'Ortenzio, F. and Ribera d'Alcalà, M.: On the trophic regimes of the Mediterranean Sea: a satellite analysis, Biogeosciences, 6, 139–148,, 2009. 

Fantozzi, L., Manca, G., Ammoscato, I., Pirrone, N., and Sprovieri, F.: The cycling and sea-air exchange of mercury in the waters of the Eastern Mediterranean during the 2010 MED-OCEANOR cruise campaign, Sci. Total Environ., 448, 151–162,, 2013. 

Ferrara, R., Ceccarini, C., Lanzillotta, E., Gårdfeldt, K., Sommar, J., Horvat, M., Logar, M., Fajon, V., and Kotnik, J.: Profiles of dissolved gaseous mercury concentration in the Mediterranean seawater, Atmos. Environ., 37, 85–92,, 2003. 

Gårdfeldt, K., Sommar, J., Ferrara, R., Ceccarini, C., Lanzillotta, E., Munthe, J., Wängberg, I., Lindqvist, O., Pirrone, N., Sprovieri, F., Pesenti, E., and Strömberg, D.: Evasion of mercury from coastal and open waters of the Atlantic Ocean and the Mediterranean Sea, Atmos. Environ., 37, 73–84,, 2003. 

Gencarelli, C. N., De Simone, F., Hedgecock, I. M., Sprovieri, F., and Pirrone, N.: Development and application of a regional-scale atmospheric mercury model based on WRF/Chem: A Mediterranean area investigation, Environ. Sci. Pollut. Res., 21, 4095–4109,, 2014. 

Gosnell, K. J. and Mason, R. P.: Mercury and methylmercury incidence and bioaccumulation in plankton from the central Pacific Ocean, Mar. Chem., 177, 772–780,, 2015. 

Gosnell, K. J., Balcom, P. H., Tobias, C. R., Gilhooly, W. P., and Mason, R. P.: Spatial and temporal trophic transfer dynamics of mercury and methylmercury into zooplankton and phytoplankton of Long Island Sound, Limnol. Oceanogr., 62, 1122–1138,, 2017. 

Hammerschmidt, C. R., Finiguerra, M. B., Weller, R. L., and Fitzgerald, W. F.: Methylmercury Accumulation in Plankton on the Continental Margin of the Northwest Atlantic Ocean, Environ. Sci. Technol., 47, 3671–3677, 2013. 

Harding, G., Dalziel, J., and Vass, P.: Bioaccumulation of methylmercury within the marine food web of the outer Bay of Fundy, Gulf of Maine, PLoS One, 13, 1–30,, 2018. 

Harmelin-Vivien, M., Cossa, D., Crochet, S., Bănaru, D., Letourneur, Y., and Mellon-Duval, C.: Difference of mercury bioaccumulation in red mullets from the north-western Mediterranean and Black seas, Mar. Pollut. Bull., 58, 679–685,, 2009. 

Heimbürger, L.-E., Cossa, D., Marty, J.-C., Migon, C., Averty, B., Dufour, A., and Ras, J.: Methylmercury distributions in relation to the presence of nano- and picophytoplankton in an oceanic water column (Ligurian Sea, North-western Mediterranean), Geochim. Cosmochim. Ac., 74, 5549–5559,, 2010. 

Hines, M. E., Horvat, M., Faganeli, J., Bonzongo, J. C. J., Barkay, T., Major, E. B., Scott, K. J., Bailey, E. A., Warwick, J. J., and Lyons, W. B.: Mercury biogeochemistry in the Idrija River, Slovenia, from above the mine into the Gulf of Trieste, Environ. Res., 83, 129–139,, 2000. 

Horvat, M., Covelli, S., Faganeli, J., Logar, M., Mandić, V., Rajar, R., Širca, A., and Žagar, D.: Mercury in contaminated coastal environments; a case study: The Gulf of Trieste, Sci. Total Environ., 237/238, 43–56,, 1999. 

Horvat, M., Kotnik, J., Logar, M., Fajon, V., Zvonarić, T., and Pirrone, N.: Speciation of mercury in surface and deep-sea waters in the Mediterranean Sea, Atmos. Environ., 37, 93–108,, 2003. 

Horvat, M., Degenek, N., Lipej, L., Snoj Tratnik, J., and Faganeli, J.: Trophic transfer and accumulation of mercury in ray species in coastal waters affected by historic mercury mining (Gulf of Trieste, northern Adriatic Sea), Environ. Sci. Pollut. Res., 21, 4163–4176,, 2014. 

Jiskra, M., Heimbürger-Boavida, L.-E., Desgranges, M.-M., Petrova, M. V., Dufour, A., Ferreira-Araujo, B., Masbou, J., Chmeleff, J., Thyssen, M., Point, D., and Sonke, J. E.: Mercury stable isotopes constrain atmospheric sources to the ocean, Nature, 597, 678–682,, 2021. 

Jonsson, S., Mazrui, N. M., and Mason, R. P.: Dimethylmercury Formation Mediated by Inorganic and Organic Reduced Sulfur Surfaces, Sci. Rep., 6, 27958,, 2016. 

Jørgensen, S. E., Friis, M. B., Henriksen, J., Jørgensen, L. A., and Mejer, H. F.: Handbook of Envirnmental Data and Ecological Parameters, edited by: Jørgensen, S. E., Pergamon Press, International Society for Ecological Modelling, ebook ISBN: 9781483188782, 1979. 

Kotnik, J., Horvat, M., Tessier, E., Ogrinc, N., Monperrus, M., Amouroux, D., Fajon, V., Gibičar, D., Žižek, S., Sprovieri, F., and Pirrone, N.: Mercury speciation in surface and deep waters of the Mediterranean Sea, Mar. Chem., 107, 13–30,, 2007. 

Kotnik, J., Horvat, M., Ogrinc, N., Fajon, V., Žagar, D., Cossa, D., Sprovieri, F., and Pirrone, N.: Mercury speciation in the Adriatic Sea, Mar. Pollut. Bull., 96, 136–148,, 2015. 

Lamborg, C. H., Hammerschmidt, C. R., Bowman, K. L., Swarr, G. J., Munson, K. M., Ohnemus, D. C., Lam, P. J., Heimbürger, L.-E., Rijkenberg, M. J. A., and Saito, M. A.: A global ocean inventory of anthropogenic mercury based on water column measurements, Nature, 512, 65–68,, 2014. 

Lamborg, C. H., Hammerschmidt, C. R., and Bowman, K. L.: An examination of the role of particles in oceanic mercury cycling, Philos. T. R. Soc. A, 374, 20150297,, 2016. 

Lanzillotta, E., Ceccarini, C., and Ferrara, R.: Photo-induced formation of dissolved gaseous mercury in coastal and offshore seawater of the Mediterranean basin, Sci. Total Environ., 300, 179–187,, 2002. 

Lazzari, P., Teruzzi, A., Salon, S., Campagna, S., Calonaci, C., Colella, S., Tonani, M., and Crise, A.: Pre-operational short-term forecasts for Mediterranean Sea biogeochemistry, Ocean Sci., 6, 25–39,, 2010. 

Lazzari, P., Solidoro, C., Ibello, V., Salon, S., Teruzzi, A., Béranger, K., Colella, S., and Crise, A.: Seasonal and inter-annual variability of plankton chlorophyll and primary production in the Mediterranean Sea: A modelling approach, Biogeosciences, 9, 217–233,, 2012. 

Lazzari, P., Mattia, G., Solidoro, C., Salon, S., Crise, A., Zavatarelli, M., Oddo, P., and Vichi, M.: The impacts of climate change and environmental management policies on the trophic regimes in the Mediterranean Sea: Scenario analyses, J. Mar. Syst., 135, 137–149,, 2014. 

Lazzari, P., Solidoro, C., Salon, S., and Bolzon, G.: Spatial variability of phosphate and nitrate in the Mediterranean Sea: A modelling approach, Deep-Sea Res. Pt. I, 108, 39–52,, 2016. 

Lazzari, P., Álvarez, E., Terzić, E., Cossarini, G., Chernov, I., D'ortenzio, F., and Organelli, E.: Cdom spatiotemporal variability in the mediterranean sea: A modelling study, J. Mar. Sci. Eng., 9, 1–18,, 2021. 

Lee, C. and Fisher, N. S.: Methylmercury uptake by diverse marine phytoplankton, Limnol. Oceanogr., 61, 1626–1639,, 2016. 

Lee, C., Lee, C., and Fisher, N. S.: Bioaccumulation of methylmercury in a marine diatom and the influence of dissolved organic matter, Mar. Chem., 1–10,, 2017. 

Lehnherr, I., St. Louis, V. L., Hintelmann, H., and Kirk, J. L.: Methylation of inorganic mercury in polar marine waters, Nat. Geosci., 4, 298–302,, 2011. 

Liu, M., Zhang, Q., Maavara, T., Liu, S., Wang, X., and Raymond, P. A.: Rivers as the largest source of mercury to coastal oceans worldwide, Nat. Geosci., 14, 672–677,, 2021. 

Luo, H., Cheng, Q., and Pan, X.: Photochemical behaviors of mercury (Hg) species in aquatic systems: A systematic review on reaction process, mechanism, and influencing factor, Sci. Total Environ., 720, 137540,, 2020. 

Mahlmann, D. M., Jahnke, J., and Loosen, P.: Rapid determination of the dry weight of single, living cyanobacterial cells using the Mach-Zehnder double-beam interference microscope, Eur. J. Phycol., 43, 355–364,, 2008. 

Margirier, F., Testor, P., Heslop, E., Mallil, K., Bosse, A., Houpert, L., Mortier, L., Bouin, M. N., Coppola, L., D'Ortenzio, F., Durrieu de Madron, X., Mourre, B., Prieur, L., Raimbault, P., and Taillandier, V.: Abrupt warming and salinification of intermediate waters interplays with decline of deep convection in the Northwestern Mediterranean Sea, Sci. Rep., 10, 1–11,, 2020. 

Mason, R. P., Choi, A. L., Fitzgerald, W. F., Hammerschmidt, C. R., Lamborg, C. H., Soerensen, A. L., and Sunderland, E. M.: Mercury biogeochemical cycling in the ocean and policy implications, Environ. Res., 119, 101–117,, 2012. 

McGeer, J. C., Brix, K. V., Skeaff, J. M., Deforest, D. K., Brigham, S. I., Adams, W. J., and Green, A.: Inverse relationship between bioconcentration factor and exposure concentration for metals: Implications for hazard assessment of metals in the aquatic environment, Environ. Toxicol. Chem., 22, 1017–1037,<1017:IRBBFA>2.0.CO;2, 2003. 

Melaku Canu, D., Rosati, G., Solidoro, C., Heimbürger, L.-E., and Acquavita, A.: A comprehensive assessment of the mercury budget in the Marano–Grado Lagoon (Adriatic Sea) using a combined observational modeling approach, Mar. Chem., 177, 742–752,, 2015. 

Mihanović, H., Vilibić, I., Šepić, J., Matić, F., Ljubešić, Z., Mauri, E., Gerin, R., Notarstefano, G., and Poulain, P. M.: Observation, Preconditioning and Recurrence of Exceptionally High Salinities in the Adriatic Sea, Front. Mar. Sci., 8, 1–22,, 2021. 

Monperrus, M., Tessier, E., Amouroux, D., Leynaert, A., Huonnic, P., and Donard, O. F. X.: Mercury methylation, demethylation and reduction rates in coastal and marine surface waters of the Mediterranean Sea, Mar. Chem., 107, 49–63,, 2007. 

Motta, L. C., Blum, J. D., Johnson, M. W., Umhau, B. P., Popp, B. N., Washburn, S. J., Drazen, J. C., Benitez-Nelson, C. R., Hannides, C. C. S., Close, H. G., and Lamborg, C. H.: Mercury Cycling in the North Pacific Subtropical Gyre as Revealed by Mercury Stable Isotope Ratios, Global Biogeochem. Cy., 33, 777–794,, 2019. 

Munson, K. M., Lamborg, C., Swarr, G. J., and Saito, M. A.: Mercury species concentrations and fluxes in the Central Tropical Pacific Ocean, Global Biogeochem. Cy., 29, 656–676, 2015. 

Munson, K. M., Lamborg, C. H., Boiteau, R. M., and Saito, M. A.: Dynamic mercury methylation and demethylation in oligotrophic marine water, Biogeosciences, 15, 6451–6460,, 2018. 

Nerentorp Mastromonaco, M. G., Gårdfeldt, K., and Wängberg, I.: Seasonal and spatial evasion of mercury from the western Mediterranean Sea, Mar. Chem., 193, 34–43,, 2017. 

Nightingale, P. D., Malin, G., Law, C. S., Watson, A. J., Liss, P. S., Liddicoat, M. I., Boutin, J., and Upstill-Goddard, R. C.: In situ evaluation of air-sea gas exchange parameterizations using novel conservative and volatile tracers, Global Biogeochem. Cy., 14, 373–387,, 2000. 

Ortiz, V. L., Mason, R. P., and Evan Ward, J.: An examination of the factors influencing mercury and methylmercury particulate distributions, methylation and demethylation rates in laboratory-generated marine snow, Mar. Chem., 177, 753–762,, 2015. 

Palanques, A., Guillén, J., Puig, P., and Grimalt, J. O.: Effects of flushing flows on the transport of mercury-polluted particulate matter from the Flix Reservoir to the Ebro Estuary, J. Environ. Manage., 260, 1–11,, 2020. 

Panagos, P., Jiskra, M., Borrelli, P., Liakos, L., and Ballabio, C.: Mercury in European topsoils: Anthropogenic sources, stocks and fluxes, Environ. Res., 201, 111556,, 2021. 

Pinardi, N. and Masetti, E.: Variability of the large scale general circulation of the Mediterranean Sea from observations and modelling: A review, Palaeogeogr. Palaeocl., 158, 153–173,, 2000. 

Pinardi, N., Zavatarelli, M., Adani, M., Coppini, G., Fratianni, C., Oddo, P., Simoncelli, S., Tonani, M., Lyubartsev, V., Dobricic, S., and Bonaduce, A.: Mediterranean Sea large-scale low-frequency ocean variability and water mass formation rates from 1987 to 2007: A retrospective analysis, Prog. Oceanogr., 132, 318–332,, 2015. 

Pinardi, N., Cessi, P., Borile, F., and Wolfe, C. L. P.: The Mediterranean sea overturning circulation, J. Phys. Oceanogr., 49, 1699–1721,, 2019. 

Ramondenc, S., Lombard, F., Santinelli, C., Stemmann, L., Gorsky, G., and Guidi, L.: An initial carbon export assessment in the Mediterranean Sea based on drifting sediment traps and the Underwater Vision Profiler data sets, Deep-Sea Res. Pt. I, 117, 107–119,, 2016. 

Rimondi, V., Costagliola, P., Lattanzi, P., Morelli, G., Cara, G., Cencetti, C., Fagotti, C., Fredduzzi, A., Marchetti, G., Sconocchia, A., and Torricelli, S.: A 200 km-long mercury contamination of the Paglia and Tiber floodplain: Monitoring results and implications for environmental management, Environ. Pollut., 255, 113191,, 2019. 

Rosati, G., Heimbürger, L. E., Melaku Canu, D., Lagane, C., Laffont, L., Rijkenberg, M. J. A., Gerringa, L. J. A., Solidoro, C., Gencarelli, C. N., Hedgecock, I. M., De Baar, H. J. W., and Sonke, J. E.: Mercury in the Black Sea: new insights from measurements and numerical modeling, Global Biogeochem. Cy., 32, 1–22,, 2018. 

Rosati, G., Solidoro, C., and Canu, D.: Mercury dynamics in a changing coastal area over industrial and postindustrial phases: Lessons from the Venice Lagoon, Sci. Total Environ., 743, 1–15,, 2020. 

Rosati, G., Canu, D., Lazzari, P., and Solidoro, C.: OGSTM-BFM-Hg model code, Zenodo [code],, 2022a. 

Rosati, G., Laurent, C., Aveytua-Alcazar, L., Solidoro, C., and Canu, D.: Investigating Hg cycling and bioaccumulation in a temperate macro-tidal lagoon impacted by legacy pollution, the Venice Lagoon case study, in: Proceedings of the 15th International Conference on Mercury as a Global Pollutant, ICMGP, Online conference, 2022b. 

Salon, S., Cossarini, G., Bolzon, G., Feudale, L., Lazzari, P., Teruzzi, A., Solidoro, C., and Crise, A.: Novel metrics based on biogeochemical argo data to improve the model uncertainty evaluation of the cmems mediterranean marine ecosystem forecasts, Ocean Sci., 15, 997–1022,, 2019. 

Schaefer, J. K. and Morel, F. M. M.: High methylation rates of mercury bound to cysteine by Geobacter sulfurreducens, Nat. Geosci., 2, 123–126, 2009. 

Schartup, A. T., Ndu, U., Balcom, P. H., Mason, R. P., and Sunderland, E. M.: Contrasting effects of marine and terrestrially derived dissolved organic matter on mercury speciation and bioavailability in seawater, Environ. Sci. Technol., 49, 5965–5972,, 2015. 

Schartup, A. T., Qureshi, A., Dassuncao, C., Thackray, C. P., Harding, G., and Sunderland, E. M.: A Model for Methylmercury Uptake and Trophic Transfer by Marine Plankton, Environ. Sci. Technol., 52, 654–662,, 2018. 

Schlitzer, R.: Ocean Data View, (last access: 26 July 2022), 2014. 

Sonke, J. E., Heimbürger, L. E., and Dommergue, A.: Mercury biogeochemistry: Paradigm shifts, outstanding issues and research needs, Comptes Rendus-Geosci., 345, 213–224,, 2013. 

Sunderland, E. M., Krabbenhoft, D. P., Moreau, J. W., Strode, S. A., and Landing, W. M.: Mercury sources, distribution, and bioavailability in the North Pacific Ocean: Insights from data and models, Global Biogeochem. Cy., 23, 1–14,, 2009. 

Terzić, E., Lazzari, P., Organelli, E., Solidoro, C., Salon, S., D'Ortenzio, F., and Conan, P.: Merging bio-optical data from Biogeochemical-Argo floats and models in marine biogeochemistry, Biogeosciences, 16, 2527–2542,, 2019. 

Tseng, C. M., Ang, S. J., Chen, Y. S., Shiao, J. C., Lamborg, C. H., He, X., and Reinfelder, J. R.: Bluefin tuna reveal global patterns of mercury pollution and bioavailability in the world's oceans, P. Natl. Acad. Sci. USA, 118, 1–6,, 2021. 

UNEP: Global Mercury Assessment 2018. UN Environment Programme, Chemicals and Health Branch Geneva, Switzerland, (last access: 25 July 2022), 2019. 

Vichi, M., Lovato, T., Lazzari, P., Cossarini, G., Gutierrez Mlot, E., Mattia, G., Masina, S., McKiver, W., Pinardi, N., Solidoro, C., Tedesco, L., and Zavatelli, M.: The Biogeochemical Flux Model (BFM): Equation Description and User Manual, BFM version 5.2. BFM Report series N. 1, Release 1.2, June 2020, Bologna, Italy,, 1–104, 2015.  

Vignati, D. A. L., Burdino, E., Congiu, A. M., Cicala, F., Pardos, M., Nieddu, G. F., and Ugazio, G.: Quality evaluation of sediments from 24 tributaries of the Po River, Italy, Water. Air. Soil Pollut., 190, 129–141,, 2008. 

Walsh, J. P. and Nittrouer, C. A.: Understanding fine-grained river-sediment dispersal on continental margins, Mar. Geol., 263, 34–45,, 2009. 

Wang, K., Munson, K. M., Armstrong, D. A., Macdonald, R. W., and Wang, F.: Determining seawater mercury methylation and demethylation rates by the seawater incubation approach: A critique, Mar. Chem., 219, 103753,, 2020. 

Wu, P., Zakem, E. J., Dutkiewicz, S., and Zhang, Y.: Biomagnification of Methylmercury in a Marine Plankton Ecosystem, Environ. Sci. Technol., 54, 5446–5455,, 2020. 

Wu, P., Dutkiewicz, S., Monier, E., and Zhang, Y.: Bottom-Heavy Trophic Pyramids Impair Methylmercury Biomagnification in the Marine Plankton Ecosystems, Environ. Sci. Technol., 55, 15476–15483,, 2021. 

Žagar, D., Petkovšek, G., Rajar, R., Sirnik, N., Horvat, M., Voudouri, A., Kallos, G., and Četina, M.: Modelling of mercury transport and transformations in the water compartment of the Mediterranean Sea, Mar. Chem., 107, 64–88,, 2007. 

Žagar, D., Sirnik, N., Četina, M., Horvat, M., Kotnik, J., Ogrinc, N., Hedgecock, I. M., Cinnirella, S., Pirrone, Nicola, De Simone, F., and Gencarelli, C. N.: Mercury in the Mediterranean, Part 2: Processes and mass balance, Environ. Sci. Pollut. Res., 21, 4081–4094,, 2014. 

Zhang, C., Dang, H., Azam, F., Benner, R., Legendre, L., Passow, U., Polimene, L., Robinson, C., Suttle, C. A., and Jiao, N.: Evolving paradigms in biological carbon cycling in the ocean, Natl. Sci. Rev., 5, 481–499,, 2018. 

Zhang, L., Wu, S., Zhao, L., Lu, X., Pierce, E. M., and Gu, B.: Mercury Sorption and Desorption on Organo-Mineral Particulates as a Source for Microbial Methylation, Environ. Sci. Technol., 53, 2426–2433,, 2019. 

Zhang, Y., Jaeglé, L., and Thompson, L.: Natural biogeochemical cycle of mercury in a global three-dimensional ocean tracer model, Global Biogeochem. Cy., 28, 553–570,, 2014a. 

Zhang, Y., Jaeglé, L., Thompson, L. A., and Streets, D. G.: Six centuries of changing oceanic mercury, Global Biogeochem. Cy., 28, 1251–1261,, 2014b. 

Zhang, Y., Soerensen, A. L., Schartup, A. T., and Sunderland, E. M.: A global model for methylmercury formation and uptake at the base of marine food webs, Global Biogeochem. Cy., 34, 1–21,, 2020. 

Short summary
Methylmercury (MeHg) is produced and bioaccumulated in marine food webs, posing concerns for human exposure through seafood consumption. We modeled and analyzed the fate of MeHg in the lower food web of the Mediterranean Sea. The modeled spatial–temporal distribution of plankton bioaccumulation differs from the distribution of MeHg in surface water. We also show that MeHg exposure concentrations in temperate waters can be lowered by winter convection, which is declining due to climate change.
Final-revised paper