Dissolved oxygen budget in the Levantine Sea: a coupled physical-biogeochemical modelling approach

The Levantine Basin is an ultra-oligotrophic region and the formation site of Levantine Intermediate Water. A high-resolution 3D coupled hydrodynamic-biogeochemical model (SYMPHONIE-Eco3MS) was used to investigate the seasonal and interannual variability of dissolved oxygen (O₂) in the Levantine Basin and to estimate its basin-wide budget over the period 2013–2020. The model results show a pronounced seasonal cycle of air–sea exchanges. During winter, cooling and vertical mixing induce an undersaturation in oxygen of the surface layer by up to 2 % across the entire basin, leading to atmospheric oxygen absorption. In contrast, during the stratified period, primary production and warming induce a slight oversaturation and subsequent oxygen release to the atmosphere. The annual budget over the 7-year period shows that the basin acts as a net sink for atmospheric oxygen. The oxygen budget analyses further indicate that the surface layer (0–150 m) acts as a source of dissolved oxygen for intermediate depths through winter vertical export, whose amplitude is significantly governed by the magnitude of heat fluxes. At the basin and annual scale, we estimate a net lateral oxygen input into the basin from the Ionian Sea and a net export towards the Aegean Sea, with this lateral export at both surface and intermediate layers enhanced when winter heat loss is intense. Biogeochemically, the Levantine Basin alternates between autotrophic and heterotrophic states on an annual basis, depending on the intensity of winter surface heat loss. Spatially, the Rhodes Gyre, a quasi-permanent cyclonic structure and major site of intermediate water formation, emerges as a significant oxygen pump in winter, with annual uptake rates twice as high as the rest of the Levantine Basin, and shows enhanced biological production during the productive season, contributing to 41 % of the net annual oxygen production in the surface layer in the basin. This study highlights the need for further modeling studies on pluri-annual and multi-decadal scales to explore interannual variability and evolution of the annual oxygen budget across the entire Eastern Basin, particularly in the context of climate change.

1 Introduction

Dissolved oxygen (O₂) is essential for marine life, supporting respiration of living organisms and the oxidation of organic matter, thereby regulating nutrient cycling and organic matter remineralization, and influencing the biogeochemical cycles of important elements in the ocean (Breitburg et al., 2018; Gruber, 2011; Morée et al., 2023). The ocean's oxygen inventory is primarily controlled by its production through photosynthesis, and its consumption through remineralization, as well as by temperature and salinity-dependent oxygen solubility, air–sea exchange and the mixing and advective fluxes influencing the ventilation of water masses (Sanders et al. 2026; Helm et al., 2011). Since 1960, the total oxygen inventory has decreased by 2 % in the Global Ocean (Schmidtko et al., 2017), a decline primarily attributed to warming-induced reductions in oxygen solubility and to enhanced upper-ocean stratification, which limits vertical ventilation (Helm et al., 2011; Schmidtko et al., 2017; Breitburg et al., 2018; Stramma and Schmidtko, 2021). However, oxygen changes present large regional and temporal heterogeneity (Schmidtko et al., 2017; Levin, 2018; Feucher et al., 2022; Kolodziejczyk et al., 2024; Wu et al., 2025), making long-term trends difficult to detect in the upper ocean. Identifying the relative contribution of physical and biogeochemical drivers is therefore essential to better understand regional oxygen dynamics.

The Levantine Basin, in the south-easternmost Mediterranean Sea (Fig. 1) is an ultra-oligotrophic region characterized by exceptionally low primary productivity (Kress and Herut, 2001) and is particularly sensitive to changes in ventilation and biogeochemical processes. As the area of the formation of the Levantine Intermediate Water (LIW) which subsequently supplies the Eastern Intermediate Water (EIW), propagating throughout the entire Mediterranean Sea at intermediate depths (Brasseur et al., 1996), the Levantine Basin plays a crucial role in basin-scale ventilation (Kress et al., 2003). Its vertical structure of dissolved oxygen reflects the complex interplay of physical, biological, and chemical processes occurring at several temporal and spatial scales. The surface and intermediate waters characterized by exchanges with the atmosphere are well oxygenated, with a pronounced ventilation in the Rhodes Gyre, a permanent cyclone in the northwest of the basin, which has traditionally been identified as the major area of LIW formation (Lascaratos et al., 1999; Lascaratos and Nittis, 1998; Sur et al., 1993). The upper layer (0–150 m) exhibits seasonal variability, with maximum oxygen values located in the mixed layer in winter, while during the stratification period, a subsurface oxygen maximum layer develops near 80 m depth, mostly attributed to both physical trapping of oxygen in Atlantic Water and biological production (Kress and Herut, 2001; Di Biagio et al., 2022). Beyond the seasonal signal associated with local processes, the properties of the upper layer are also influenced by the general eastern Mediterranean circulation and in particular the Adriatic-Ionian Bimodal Oscillation System (BiOS) (Gacic et al., 2010, 2011; Velaoras et al. 2014; Menna et al., 2022), characterized by two alternating circulation regimes: during the anti-cyclonic phase of the Northern Ionian Gyre (NIG), Atlantic Water (AW) flowing across the Sicily Channel is preferentially directed northwards, toward the Adriatic Sea, while during its cyclonic phase, AW, directed eastward, mostly supplies the Levantine Basin. The BiOS process presents quasi-decadal variability and has been proposed as a driver of changes in the thermohaline and biogeochemistry in the Adriatic and Levantine seas (Civitarese et al., 2010; Velaoras et al. 2014; Ozer et al., 2017, 2022; Di Biagio et al., 2023; Civitarese et al., 2023). Below the euphotic layer, oxygen concentrations decline with an Oxygen Minimum Layer (OML) located below the Intermediate Water between 600 and 1200 m (Cardin et al., 2015; Mavropoulou et al., 2020) and characterized by concentrations of 170/180 µmol kg⁻¹ (Tanhua et al., 2013). Deep water masses of the sub-basin are sensitive to variations in deep-water formation and circulation taking place in the Eastern Mediterranean. After the Eastern Mediterranean Transient (EMT) in the early 1990's when the deep water formation area shifted from the Adriatic to the Aegean Sea, an increase of oxygen was documented in the deeper layers in response to the inflow below 2500 m of more oxygenated waters originating from the Aegean Sea (Lascaratos et al., 1999; Mavropoulou et al., 2020). This was accompanied with the upward displacement of the older Adriatic-origin deep waters and of the OML. More recent observations, however, indicate a deoxygenation trend since 2008, attributed to weakened deep water formation and reduced ventilation rates as a progressive return to pre-EMT characteristics with an homogenized deep layer (Sisma-Ventura et al., 2021).

Figure 1 Model domain and bathymetry (m, background) in the Eastern Mediterranean. The arrows represent the simulated surface currents averaged over the study 7-year period (2013–2020), black thick lines delimit the basin for the budget calculation. Red, yellow, and green dots indicate BGC-Argo floats trajectories, PERLE-1, and PERLE-2 cruise stations, respectively, and white crosses CARIMED cruise stations, over the period from 2013 to 2021. Blue dots in the insert represent the river mouths.

To date, with limited spatial and temporal observations in the area, the variability of oxygen inventory in the Levantine Basin remains poorly understood, and there is no proposed comprehensive budget quantification for the entire region. In the framework of the PERLE (Pelagic Ecosystem Response to Deep Water Formation in the Levant Experiment, Conan and Durrieu De Madron, 2019) project, the present work aims at quantifying the seasonal and interannual variations in the oxygen inventory of the Levantine surface and intermediate water masses, detailing the contribution of air–sea oxygen fluxes, biogeochemical and physical fluxes. This analysis is based on 3D coupled hydrodynamic-biogeochemical model outputs covering a period of 7 years, from 2013 to 2020. Following on the budget approach developed by Ulses et al. (2021) for the north-western Mediterranean Sea, we investigate the ultra-oligotrophic Levantine Basin and provide a basin-scale quantification of its dissolved oxygen budget, highlighting the role of transport processes and permanent circulation features such as the Rhodes Gyre.

After the introduction (Sect. 1), this paper is organized as follows. Section 2 describes the coupled hydrodynamic-biogeochemical model implemented in the Levantine Basin and an assessment of the model results using in situ observations. Section 3 investigates the seasonal and interannual dynamics of oxygen in the surface and intermediate layers for the Levantine Basin, estimates an annual budget of oxygen, and finally describes its spatial variability. This section is followed by a discussion of the results and a conclusion in Sects. 4 and 5, respectively.

2 Material and Method

2.1 Modeling

2.1.1 The coupled hydrodynamic-biogeochemical model

This study is based on a Mediterranean configuration of the ocean circulation model SYMPHONIE (Marsaleix et al. 2006, 2008), forcing offline the biogeochemical model Eco3M-S (Ulses et al., 2016, 2023). The horizontal resolution of the model grid varies from 2.3 to 4.5 km, with a refined resolution of 1.3 km in the Gibraltar Strait. The vertical grid has 60 vertical vanishing quasi sigma levels. More details can be found in Estournel et al. (2021).

Eco3M-S is a multi-nutrient and multi-plankton functional type model, representing the dynamics of the pelagic plankton ecosystem and the cycles of carbon, nitrogen, phosphorus, silicon, and oxygen (Auger et al., 2011; Ulses et al., 2023), with 37 state variables. The rate of change of dissolved oxygen concentration due to biogeochemistry in the water column is calculated based on the following equation:

$$\begin{array}{}\text{(1)}& \begin{array}{rl}& {\displaystyle \frac{\mathrm{dDOx}}{\mathrm{d}t}}=\sum _{i=\mathrm{1}}^{\mathrm{3}}\left({\mathrm{GPP}}_i\right.\left.-{\mathrm{RespPhy}}_i\right){\mathit{\gamma }}_{\mathrm{C}/\mathrm{DOx}}\\ & -\sum _{i=\mathrm{1}}^{\mathrm{3}}\left({\mathrm{RespZoo}}_i\right){\mathit{\gamma }}_{\mathrm{C}/\mathrm{DOx}}-\mathrm{RespBac}{\mathit{\gamma }}_{\mathrm{C}/\mathrm{DOx}}\\ & +({\mathrm{UptPhy}}_{i,{\mathrm{NO}}_{\mathrm{3}}}-\mathrm{Nitrif}){\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}/\mathrm{DOx}}\end{array}\end{array}$$

The dissolved oxygen concentration is represented by the term DOx. GPP_i and RespPhy_i are gross primary production and respiration, respectively, for phytoplankton group i. RespZoo_i and RespBac are respiration of zooplankton group i and of bacteria, UptPhy${}_{i,{\mathrm{NO}}_{\mathrm{3}}}$, and Nitrif uptake of nitrate by phytoplankton class i, and nitrification, respectively. ${\mathit{\gamma }}_{\mathrm{C}/\mathrm{DOx}}$ and ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}/\mathrm{DOx}}$ , equal to 1 and 2, respectively, are the moles of DOx used per mole of C in respiration and needed to oxidize one mole of ammonium in nitrification as described in Grégoire et al. (2008). The flux of dissolved oxygen at the air–sea interface is calculated using the Garcia and Gordon (1992) equation for the solubility, and the parametrization of Wanninkhof and McGillis (1999) for the transfer velocity, following the study of Ulses et al. (2021) in the northwestern Mediterranean deep convection area.

2.1.2 Initialisation and boundary conditions

The implementation of the coupled model was described in detail in Estournel et al. (2021) and Habib et al. (2023). The period simulated by the hydrodynamic model runs from May 2011 to May 2021. Atmospheric forcings for both hydrodynamic and biogeochemical models were provided by the HRES (atmospheric model High RESolution forecast) product of ECMWF model with a horizontal resolution of $\mathrm{1}/\mathrm{8}$° using hourly fields (wind, air temperature and humidity, pressure, solar and downward longwave radiation, and precipitation). The biogeochemical model Eco3M-S was forced by daily fields of temperature, salinity, current, and vertical diffusivity from the SYMPHONIE model. It covers the period between August 2011 till March 2021. The first two years (July 2011–December 2013) of the biogeochemical simulation were dedicated to model spin-up to ensure biogeochemical stability and were not considered in the analysis, while the December 2013–December 2020 period was used for the budget analysis. This period was selected based on the availability of consistent physical forcing and the density of in situ observations for model initialization and validation. The biogeochemical model was initialized using climatological fields of in situ nutrient and dissolved oxygen concentrations from the CARIMED (CARbon in the MEDiterranean Sea, Álvarez et al., 2019) database and Biogeochemical-Argo (BGC-Argo) float data over the 2011–2012 summer periods when data were available, in 10 sub-regions. At the mouths of 142 rivers taking into account (Fig. 1), concentrations of nutrients were imposed by sub-basin using the dataset of Ludwig et al. (2010). Dissolved oxygen at river mouths was set at saturation values. In the Atlantic Ocean, nutrients were prescribed using monthly profiles from the World Ocean Atlas 2009 climatology at 5.5° W. In the Marmara Sea, to represent a two-layer flow regime, we imposed a daily relaxation towards a nutrient concentration of 0.24 and 1.03 mmol N m⁻³ and a phosphate concentration of 0.06 and 0.05 mmol P m⁻³ for depths above and below 15 m, respectively, based on the observations near the Dardanelles Strait from (Tugrul et al., 2002).

2.1.3 Budget calculation

For spatial mean and budget calculation, we defined an area (delimited by the black lines in Fig. 1) covering 540 000 m², with the boundary with the Ionian Sea linking the southwestern Cretan coast to the Libyan coast. The water column was divided into three layers based on the thermohaline structure represented by the physical model (Estournel et al., 2021) and the associated dominant biogeochemical processes: the surface layer defined as the photic layer covering the surface to 150 m depth where photosynthesis takes place, the underlying intermediate layer from 150 to 400 m where LIW flows, and the deep layer below 400 m. In this study, we will be focusing on the first two layers where changes, in particular related to LIW formation, occur generally more rapidly. The biogeochemical term of the oxygen budget is the sum of oxygen production due to gross primary production and nitrate uptake by phytoplankton, and of oxygen consumption through nitrification and community respiration. The physical term is divided into two components: the lateral and the vertical transports, which are both due to advection and mixing processes. The lateral transport represents the exchanges at the boundaries with the Ionian and Aegean seas. A negative lateral transport indicates a net export of oxygen from the considered layer of the Levantine Basin. The oxygen inventory, air–sea fluxes, biogeochemical fluxes, and lateral fluxes were calculated online while the vertical transport, defined as a net flux at the layer interface, was deduced from the other terms of the budget. The budget calculation is detailed in Sect. S1 in the Supplement.

2.2 Model assessment

An assessment of the hydrodynamic and biogeochemical simulations in the surface and intermediate water masses has been performed in previous studies (Estournel et al., 2021; Habib et al., 2023), suggesting its capacity to reproduce the observed general hydrology and biogeochemistry (chlorophyll, dissolved inorganic nutrients, and dissolved oxygen) in the Levantine Sea. Here, the model is further assessed in terms of the dissolved oxygen dynamics by providing supplementary comparisons with observations from BGC-Argo floats (6901528 and 6901764, http://www.coriolis.eu.org, last access: 11 February 2026), PERLE cruises (PERLE-1 and PERLE-2, https://campagnes.flotteoceanographique.fr/campagnes/18000848/fr/, last access: 5 June 2025, Conan and Durrieu De Madron, 2019), and those gathered in the CARIMED database (Álvarez et al., 2025), as well as from situ measurements of metabolic rates.

The model reproduces the observed seasonal and vertical variability (Figs. 2, S1), with strong agreement with observations throughout the upper 500 m, in particular in the timing and depth of the subsurface oxygen maximum and the oxygen minimum layer. The discrepancies found concern a slight underestimation of the subsurface oxygen maximum (RMSD ≃ 8 µmol kg⁻¹) compared to the BGC-Argo float, and an overestimation below 100 m compared to PERLE-2 observations, likely due to an overestimation of remineralization processes or vertical diffusion. The model and in situ data significantly correlate with correlation coefficient higher than 0.95 (p-value < 0.05). The RMSD values (Root Mean Square Difference) between modeled and observed surface oxygen and solubility are less than 5 µmol O₂ kg⁻¹ for both Argo floats and fall within the oxygen uncertainty interval associated with Argo float data (∼ 2–10 µmol kg⁻¹ depending on the sensor, Grégoire et al., 2021). Finally, comparisons between model results with metabolic rate measurements near the surface and within the upper layer over the May-July period (BOUM cruise (Christaki et al., 2011); THRESHOLD cruises (Regaudie-de-Gioux et al., 2009); MINOS cruise (Moutin and Raimbault, 2002) indicate that modeled GPP (gross primary production), CR (community respiration) and NCP (net community production, corresponding to GPP minus CR) fall in the observed range, generally in their upper values (Table S1). These comparisons show the model's ability to represent oxygen-related biogeochemical processes in the Levantine Basin.

Figure 2 From top to bottom: (a) trajectory of the BGC-Argo floats with deployment position (red cross) and chronology in color; Hovmöller diagrams of oxygen concentration (µmol O₂ kg⁻¹) from (b) float data and (c) model outputs for the first 500 m; (d) surface oxygen concentration in the first 10 m (µmol O₂ kg⁻¹) and (e) oxygen solubility (µmol O₂ kg⁻¹), from the float data (red) and the model (blue).

3 Results

3.1 Seasonal variability

Figure 3 presents the mean annual cycle of the modeled air–sea heat flux, wind stress, mixed layer (ML) depth, and surface temperature, spatially averaged over the Levantine Sea from December 2013 to December 2020. During fall, the decrease in air temperature leads to significant sea surface heat loss, while intensified northern winds weaken stratification, gradually deepening the mixed layer (Fig. 3a–c). The sea surface temperature drops significantly (Fig. 3d). Heat loss events persist through winter, leading the surface temperature to reach a minimum of approximately 17 °C and the mixed layer depth to gradually increase, peaking in January/February. The yearly maximum ML depth averaged spatially over the seven years is 108 ± 11 m (Table S2). In March/April, the sea surface starts gaining heat, and the surface temperature increases (Fig. 3a and d). The ML abruptly shallows but still exhibits large variations during early spring, in response to the events of continental cold winds. The frequency of intense wind events decreases in late spring/summer (Fig. 3b). Surface temperature reaches maximum values around 28 °C in August (Fig. 3d), and a thin ML settles until October. In the following, the annual cycle is divided into two successive periods based on the vertical mixing intensity. The first period is a mixing period, from October to March, and the second period is a stratified period, from April to September. The two periods were defined based on a mixed layer depth threshold of 25 m, following the criteria used by D'Ortenzio et al. (2008) and Houpert et al. (2015).

Figure 3 Annual time series of modeled (a) air–sea heat fluxes (W m⁻²), (b) wind stress (N m⁻²), (c) mixed layer depth (m), and (d) surface temperature (°C), averaged over the Levantine Sea and the period 2013–2020. In (c) and (d), the solid line corresponds to the temporal mean, the shaded area to the standard deviation. The grey shaded area represents the stratification period.

The annual cycle of the 7-year averaged oxygen fluxes and inventory variations is shown in Figs. 4 and 5 for the surface and intermediate layers, respectively. In both layers, the oxygen content increases during the mixing period and gradually decreases during the stratified period, with minimum values in November/December for the surface layer and January for the intermediate layer (Figs. 4a and 5a). Both vertical and net horizontal transports exhibit a clear seasonal variation, reaching maximum values of 50 mmol m⁻² d⁻¹ in the surface layer and 50 and 40 mmol m⁻² d⁻¹, respectively in the intermediate layer (Figs. 4f and 5e). During the mixing period, and particularly during events of strong winds and vertical mixing (Fig. 3b–c), oxygen is exported from the surface layer towards the intermediate layer (Fig. 4f), and subsequently, from the intermediate towards the deep layer (Fig. 5e). The associated mean export rates amount to 0.56 and 0.46 mol O₂ m⁻² per month, respectively (Fig. S3b, d). During the stratified period, the downward export of O₂ towards the intermediate and deeper layers is reduced by 75 % and 40 %, respectively (Fig. S3d).

Figure 4 Oxygen concentration and budget of the 0–150 m layer of the Levantine Basin averaged over the period of study. (a) Variation of the dissolved oxygen inventory (mol m⁻²) relative to initial conditions (Values are normalized to the starting time point), (b) surface oxygen saturation (orange) = (DOx − DOx sat)/DOx sat × 100 %), (c) air-sea flux (positive values correspond to downward fluxes, mmol O₂ m⁻² d⁻¹), (d) biogeochemical flux (mmol O₂ m⁻² d⁻¹), (e) sum of vertical (through the 150 m depth) and lateral (exchanges with the Ionian and Aegean Seas) transport fluxes (mmol O₂ m⁻² d⁻¹), (f) vertical (light blue) and lateral (dark blue) fluxes (mmol O₂ m⁻² d⁻¹), (g) lateral fluxes at the boundary with the Ionian (purple) and Aegean (red) Seas (mmol O₂ m⁻² d⁻¹). Horizontal transport fluxes are scaled to the area of the Levantine Basin for comparison with the other budget terms. The grey shaded area represents the stratification period.

Figure 5 Mean annual cycle of (a) variation of the dissolved oxygen inventory (mol m⁻²) relative to initial conditions (Values are normalized to the starting time point), and the different oxygen fluxes (mmol m⁻² d⁻¹): (b) biogeochemical flux, (c) total vertical and horizontal transport, (d) vertical (downward) flux (light blue) and lateral flux (dark blue), (e) the vertical fluxes at 150 and 400 m and (f) the lateral Ionian (purple) and Aegean (red) fluxes, in the intermediate layer (150–400 m) and averaged over the Levantine Basin. The grey shaded area represents the stratification period.

Horizontal oxygen transport in the surface layer is characterized by a net inflow from the Ionian Sea and an outflow towards the Aegean Sea (Figs. 4g, S2b). These exchanges are stronger during the mixing period than during the stratified period (with inflow of 1.3 vs. 0.6 mol O₂ m⁻² per month from the Ionian Sea and outflow of 0.9 vs. 0.5 mol O₂ m⁻² per month from the Aegean Sea, respectively, Fig. S3b). In the intermediate layer, horizontal exchanges are weaker overall but still display a seasonal signal. During the mixing period, a stronger net inflow from the Ionian Sea (0.2 mol O₂ m⁻² per month) and a stronger outflow towards the Aegean Sea (−0.4 O₂ mol m⁻² per month, Figs. S3d and 5f) are obtained. While the net oxygen transport in the surface layer remains directed from the Ionian and towards the Aegean across both periods, the intermediate layer exhibits net oxygen export toward both the Aegean and Ionian Seas during the stratified period, accounting for 96 % and 4 % of the total horizontal export respectively (Fig. S3d).

Model results indicate that the ecosystem in the surface layer of Levantine Basin acts as a net source of dissolved oxygen from January to August, with the highest production between February and March. In contrast, net consumption dominates from September to December (Fig. 4d). Overall, the biogeochemical O₂ flux results in a mean oxygen loss of 0.03 mol m⁻² per month during the mixing period (from October to April) and a gain of 0.06 mol m⁻² per month during the stratified period (Fig. S3a). Maximum magnitudes of biological production (> 2 mmol O₂ m⁻² d⁻¹) are located near the surface during the periods of winter mixing and the associated phytoplankton bloom, before shifting to subsurface layers later in the year (Fig. S4f). Oxygen consumption peaks in fall between 100 and 200 m depth, and during the mixing period below the mixed layer. During the stratified period, a relatively thick subsurface oxygen maximum (SOM) layer persists around 70 m depth, overlying the subsurface maximum of biological production located at around 140 m depth, close to the deep chlorophyll maximum (Fig. S4b, c, and f). This is in agreement with the findings of the modelling study of Di Biagio et al. (2022). In the intermediate layer (150–400 m), the ecosystem is characterized by a loss of oxygen throughout the year, with biogeochemical fluxes reaching values lower −5 mmol m⁻² d⁻¹ (Fig. 5b).

The air–sea oxygen flux displays a marked seasonal pattern (Fig. 4b). During the October–April mixing period, the Levantine Basin, which is undersaturated in oxygen compared to the atmosphere, acts as a sink of atmospheric oxygen. From September onward, oxygen solubility has increased (Fig. 4b) with the decrease in surface temperature (Fig. 3d) since September. In parallel, the gradual deepening of the mixed layer favors an increase in the surface oxygen concentration, through the mixing of surface O₂ poorer waters with subsurface O₂-rich waters (Fig. S4c), although surface concentrations remain below the saturation level. The air-sea flux is particularly strong in winter under strong wind conditions, and reaches maximal values around 70 mmol m⁻² d⁻¹ in early January (Fig. 4c). When averaged over the mixing period, the air–sea flux amounts to 0.50 mol O₂ m⁻² per month (Fig. S3a). At the onset of the stratified period (April–May), surface oxygen concentration increases to 230 µmol kg⁻¹, slightly exceeding saturation levels due to biological oxygen production in the surface layer (Sect. 3.2.4). As a result, the Levantine Basin becomes a source of oxygen for the atmosphere (Fig. 4c). During the rest of the stratified period, the surface oxygen concentration continues to present values higher than the oxygen solubility, leading to continuous outgassing of O₂. We estimate a mean net release of 0.26 mol O₂ m⁻² per month of oxygen to the atmosphere over the whole stratified period (Fig. S3a).

3.2 Interannual variability

The analysis of the seasonal cycle shows that oxygen fluxes exhibit their largest variability during winter, as reflected by higher standard deviations. In addition, correlations between annual oxygen processes (NCP, downward export, …) and seasonal fluxes indicate that winter is the most influential season on oxygen dynamics (Fig. S5). Therefore, the following section focuses on winter conditions to investigate the interannual variability of oxygen-related processes. Winter (December-January-February) heat loss exceeds the seasonal mean value of 152 W m⁻² for years 2014–2015, 2016–2017, 2018–2019, and 2019–2020 (Fig. 6a, Table S2). In contrast, wind stress does not show consistent interannual signals, with peak values around 0.2 N m⁻² occurring every winter (Fig. 6b). The ML depth presents interannual variability primarily associated with variations in heat loss fluxes (Fig. 6a and c). ML deeper than the 7-year mean value of 108 m are found during the winters 2014–2015, 2016–2017, and 2018–2019 (Table S2). Based on the mean winter heat flux (W-HF), and mean and maximum ML, the years were classified into two categories: as mild and cold winter years. Years with both winter heat loss and maximum ML above the seven-year means, i.e. 2014–2015, 2016–2017 and 2018–2019, are classified as cold winter years, whereas the remaining years, 2013–2014, 2015–2016, 2017–2018 and 2020–2021, are classified as mild years.

Figure 6 Time series of modeled (a) air–sea heat fluxes (W m⁻²), (b) wind stress (N m⁻²), and (c) mixed layer depth (m), averaged over the Levantine Basin. The mean winter (December-January-February) heat loss is indicated in (a). The blue shaded area represents the winter of the cold winter years.

All years display qualitatively similar seasonal oxygen cycles in the surface layer (Fig. 7). However, cold winter years (2014–2015, 2016–2017, 2018–2019) exhibit substantially larger oxygen fluxes and stronger inventory variations. During these winters, changes in oxygen inventory exceed 2 mol O₂ m⁻² (Fig. 7a). During phytoplankton blooms in these cold winters, biogeochemical fluxes surpass 20 mmol m⁻² d⁻¹ (Fig. 7d). Lateral and vertical transports of dissolved oxygen at the surface layer boundaries are also enhanced during cold winters. Lateral inflow from the Ionian Sea, lateral outflow toward the Aegean Sea, and downward export to the intermediate layer, show peak values exceeding 100, 75, and 100 mmol m⁻² d⁻¹, respectively (Fig. 7f–g). In addition, air–sea oxygen fluxes also exceed 150 mmol m⁻² d⁻¹ during cold winters (Fig. 7c).

Figure 7 As in Fig. 4, but showing the interannual variability of oxygen inventory and budget components for the 0–150 m layer of the Levantine Basin. Blue shading marks cold winters.

In the intermediate layer, cold winter years likewise display the strongest O₂ winter fluxes and the largest oxygen inventory variation compared to mild winter years (Fig. 8). Biogeochemical fluxes exhibit pronounced negative peaks exceeding −10 mmol m⁻² d⁻¹ during cold winter years (Fig. 8b). While these negative extrema are not significantly stronger than those obtained in mild years, the subsequent positive fluxes following cold winters are markedly enhanced, except during 2016–2017. Physical fluxes also are marked by larger magnitudes during cold winters (Fig. 8c–f), with enhanced lateral and vertical oxygen exchanges, as found in the upper layer. Along with the seasonal internal variation, an increasing trend in the inventory is visible from 2013–2014 to 2018–2019, followed by a decreasing trend until the end of the study period (Fig. 8a).

Figure 8 As in Fig. 5, but showing the interannual variability of oxygen inventory and budget components for the 150–400 m layer of the Levantine Basin. The blue shaded area represents the winter of the cold winter years.

Overall, oxygen flux variability is partly linked to the variability of winter heat loss (W-HL). A strong correlation is found between mean W-HL and winter downward oxygen export from the surface layer (R= 0.76, p-value <0.05, Fig. S5). Cold years also show enhanced NCP (Net Community Production: gross primary production (GPP) minus community respiration (CR)), with a significant correlation between mean W-HL and annual NCP (R = 0.91, p-value < 0.05, Fig. S5). When W-HL drops below 135 W m⁻², the trend suggests a shift of the system from autotrophic to heterotrophic conditions. In the intermediate layer, annual oxygen consumption remains relatively stable but is still significantly correlated with W-HL (R= 0.94, p-value < 0.05). Air-sea oxygen fluxes also show strong correlations with W-HL, both at the annual scale (R= 0.92, p-value < 0.05) and during winter (R= 0.93, p-value < 0.05). Although lateral fluxes intensify during cold years, only exchanges with the Aegean Sea show a significant correlation with W-HL (R= 0.74 in the surface layer, R= 0.82 in the intermediate layer). Correlations with exchanges from the Ionian Sea are weaker and not statistically significant (R= 0.69, p-value < 0.08).

3.3 Annual oxygen budget in the surface and intermediate layers

On an annual basis and averaged over the seven studied years (Fig. 9), the surface ecosystem of the Levantine Basin generally acts as a net source of oxygen, although the magnitude and even the sign of the biogeochemical contribution vary substantially from year to year. The surface layer of the basin also exhibits a net gain of oxygen through air–sea exchange. This atmospheric input results from the seasonal alternation between wintertime uptake and summertime outgassing described in Sect. 3.1. In addition to the important air–sea fluxes, the basin receives a significant lateral supply of oxygen through the advection of waters originating from the Ionian Sea into the surface layer. A fraction of the oxygen supplied to the surface is exported downward toward intermediate depths, partially consumed by biogeochemical processes, while a remaining fraction is transported toward the Aegean Sea. On average, the oxygen budget indicates that the net biogeochemical oxygen flux is one order of magnitude smaller than the contributions from transport and air–sea exchange.

Figure 9 Schematic showing the terms of the mean annual oxygen budget (in mol O₂ m⁻² yr⁻¹) for the Levantine Basin over the period from December 2013 to December 2020. The terms of the budget are estimated for the upper (surface-150 m) and intermediate (150–400 m) layers. SW: surface layer, IW: intermediate layer.

In the intermediate layer, the oxygen loss due to biogeochemical consumption is compensated by transport, at the annual scale (Fig. 9). The physical supply to the intermediate layer arises from both downward transport from the surface layer and lateral inflow from the Ionian Sea. A fraction of oxygen in intermediate waters is also exported toward the Aegean Sea (Fig. 9).

One can notice that in the surface layer, the lateral export of oxygen towards the Aegean Sea is smaller than the net input from the Ionian Sea, and on the contrary, in the intermediate layer, the lateral oxygen export toward the Aegean Sea is larger than the input from the Ionian Sea. This reflects both a small biological production in the surface layer and a large supply of oxygen from the surface layer into the intermediate waters through vertical transfer, occurring mostly during winter mixing, convection and subduction in the Levantine Sea. This budget highlights the role of the Levantine Basin as both a regional sink for intermediate waters and a transit zone for dissolved oxygen within the eastern Mediterranean circulation.

3.4 Spatial variability of oxygen fluxes in the Levantine Basin

At the annual scale, the whole Levantine Basin appears as an atmospheric sink for oxygen, except in the coastal area influenced by the Nile River (Fig. 10a). In the offshore region, the strongest oxygen uptake rates are found in the Rhodes Gyre area, covering 5 % of the surface of the Levantine Basin, but contributing 14 % of the annual atmospheric oxygen intake. Additional regions of enhanced uptake rates are located in the North, in the Antalya Bay and the Cilician Basin identified as another LIW formation area (Fach et al., 2021). The spatial distribution of the annual air-sea oxygen flux is largely controlled by the winter air–sea O₂ flux (not shown). The annual anomalies show that cold years (2014–2015, 2016–2017 and 2018–2019) are associated with enhanced atmospheric oxygen uptake over the entire basin (Table S3), and especially in the Rhodes Gyre (Fig. S8). In this main area of intermediate water formation, the vertical mixing period is characterized by the upward supply of colder and oxygen-poorer water from intermediate depths to the surface, more pronounced than in surrounding areas. The resulting negative temperature and oxygen anomaly reinforce surface undersaturation in this area, reaching maximum values varying between 2 % and 5 % (with higher values during cold winter years). In contrast to the basin-wide functioning, where oxygen is generally exported downward from surface to intermediate layers, the Rhodes Gyre displays an inverse pattern. In this area, oxygen is transported upward from intermediate depths toward the surface, and subsequently redistributed laterally. This lateral transfer is particularly strong during winter and is associated with the dispersal of LIW by subduction (Estournel et al., 2021).

Figure 10 Modeled annual air–sea oxygen flux and net community production (NCP, mol O₂ m⁻² yr⁻¹) in the surface layer (0–150 m) for the period from December 2013 to December 2020. The black line delimits the Rhodes Gyre. Positive values for the air–sea oxygen flux indicate a net flux of oxygen into the ocean (uptake), while negative values indicate a net flux from the ocean to the atmosphere. Positive NCP values correspond to net biological oxygen production and negative values to consumption.

Model results also show a marked spatial heterogeneity in the balance between GPP and CR, expressed through net community production (NCP), when averaged over the period of study in the surface layer (surface-150 m; Fig. 10b). Positive NCP values are found in a central region encompassing the cyclonic Rhodes Gyre and other cyclonic gyres (e.g., West Cyprus gyres), as well as in the coastal areas influenced by river inputs. In particular, the Rhodes Gyre contributes to 41 % of the total annual biological oxygen production in the surface layer of the entire basin. In contrast, negative values of NCP prevail in the along-slope circulation, and within the anticyclonic Mersa-Matruh Eddies and Shikmona Eddy. Annual NCP follows the same spatial pattern throughout all the years with a more pronounced production in the cyclonic gyres during cold years compared to mild winters (Fig. S7).

4 Discussion

In the present study, we used a 3D coupled physical-biogeochemical model to investigate the dynamics of oxygen in the Levantine Basin. The physical and biogeochemical parts of the coupled model were previously validated by Estournel et al. (2021) and Habib et al. (2023), respectively. Here we have further compared our results on the oxygen cycle with two types of in situ observations: high-resolution BGC-Argo data and data from research cruises. The major limitation highlighted by these comparisons is the representation of the subsurface oxygen maximum layer. The marked heterogeneity of this layer in observations, with maximum concentrations in its upper part followed by a progressive decrease of its value with depth, is not fully reproduced in the model. The increase of the maximum value during the summer period, shown in the BGC-Argo data, may reflect production or respiration processes that are underestimated or overestimated, respectively, in the model. Alternatively, this discrepancy may be explained by physical processes with a misrepresentation of the thickness of the subsurface oxygen maximum layer in the model. In that case, a finer vertical resolution at those depths or an improvement of the vertical advection scheme, avoiding possible spurious numerical mixing, as proposed by Garinet et al. (2024), could be tested in future works to reduce potential excessive vertical diffusivity.

Our budget of oxygen is subject to sources of uncertainties linked to the physical and biogeochemical models used in this study. One approach to overcome single-model uncertainties and limitations would be to adopt a multi-model approach. An alternative approach would consist in estimating oxygen budgets using observational syntheses; however, the sparse spatial coverage of in situ data currently limits the closure of basin-scale oxygen budgets based solely on observations. Finally, a complementary strategy could involve a combined approach such as that developed by Di Biagio et al. (2023) which relies on biogeochemical reanalysis corrected with Argo float data. Despite these limitations linked to the models used here, we found that the seasonal variations of oxygen solubility and concentration align with previous observational studies (Kress and Herut, 2001; Schlitzer et al., 1991). To provide a first quantitative assessment of the contribution of various processes – air–sea flux, physical dynamics, and biogeochemical processes – to oxygen budget, we chose an online and strictly closed budget approach. As shown by Trinh et al. (2024), this approach can yield substantially higher accuracy than offline calculations, especially in the quantification of lateral fluxes.

The model shows a net annual weak biological production of oxygen in the surface layer of the Levantine Basin, primarily due to the sea's oligotrophic nature, which is more pronounced in the southeastern regions of the Levantine Basin (D'Ortenzio and Ribera d'Alcalá, 2009; Lavigne et al., 2015). This oligotrophy is attributed to an anti-estuarine circulation characterized by an eastward inflow of surface nutrient-depleted waters and an outflow of intermediate nutrient-rich waters resulting from the water formation (Robinson and Golnaraghi, 1993). Interannual variability is nevertheless observed here, with the sea being heterotrophic during mild winter years (2014 and 2018), consistent with the findings by Mayot et al. (2016) using satellite ocean color data. By considering both surface and intermediate layers, the Levantine Sea appears as a net heterotrophic system in the model results. This is in line with previous studies reporting strong temporal and spatial heterogeneity in the trophic status of the oligotrophic Levantine Basin (Christaki et al., 2011; Siokou-Frangou et al., 2010). The planktonic ecosystem is largely regulated by heterotrophic processes, with higher heterotrophic/autotrophic biomass ratios typically observed in the most oligotrophic regions and during stratified periods (Christaki et al., 2002; Siokou-Frangou et al., 2002). Despite this general heterotrophic character, mesoscale physical structures promoting vertical mixing and nutrient upwelling play a key role in shaping the basin's trophic gradients, particularly within cyclonic systems (Legendre and Rassoulzadegan, 1995; Salihoglu et al., 1990). Our model results show the high contribution of the Rhodes Gyre to the annual oxygen biological production in the surface layer of the whole Levantine Basin area (around 41 %). It has been identified as the major area of LIW formation, characterising winter vertical mixing enriching the surface layer with nutrients and stimulating primary production (D'Ortenzio et al., 2021; Lavigne et al., 2013). In the intermediate layer, biogeochemical fluxes in the Gyre and the Levantine Basin exhibit little variation between mixing and stratification periods, especially during cold years, consistent with the findings of Roether and Well (2001) and Klein et al. (2003).

The model indicates that the Levantine Basin absorbs atmospheric oxygen from November to April, while releasing it during the rest of the year. This is in line with previous observational and modelling studies (Schlitzer et al., 1991; Kress and Herut, 2001; Di Biagio et al., 2022). Except for the river-influenced areas (Nile river), the whole Levantine Basin acts as an annual sink of atmospheric oxygen for all studied years, with an average uptake of 1.5 ± 0.6 mol O₂ m⁻² yr⁻¹, and higher values during cold winter years. Uptake is enhanced in intermediate water formation areas, in particular in the Rhodes Gyre, where undersaturation increases during winter due to stronger surface cooling and mixing of poorer O₂ water masses with surface waters, in agreement with what was previously observed and modeled in other water formation areas (Copin-Montégut and Bégovic, 2002; Coppola et al., 2017, 2018; Di Biagio et al., 2022; Fourrier et al., 2022; Körtzinger et al., 2004, 2008; Ulses et al., 2021; Wolf et al., 2018). The Rhodes Gyre shows a comparable winter uptake rate (20.3 ± 7.4 mol O₂ m⁻² yr⁻¹) as other water formation areas such as the Labrador Sea and Gulf of Lion (ranging between 11 and 37 mol m⁻² yr⁻¹; Copin-Montégut and Bégovic, 2002; Coppola et al., 2017, 2018; Körtzinger et al., 2008; Ulses et al., 2021; Wolf et al., 2018). As a matter of comparison, the 7-year averaged oxygen uptake estimated here for the whole Levantine Basin, characterized by relatively low solubility compared to the rest of the Mediterranean (Mavropoulou et al., 2020; Di Biagio et al., 2022), represents 64 % of the oxygen uptake by the NW Mediterranean deep convection estimated for the cold year 2012–2013 with the same coupled model (Ulses et al., 2021). These estimates are nevertheless subject to methodological uncertainties. In particular, air–sea oxygen fluxes depend on the parameterization of the gas transfer velocity, whose sensitivity to wind speed and formulation can induce uncertainties of the order of 12 %–16 %, as quantified by Ulses et al. (2021). Additional uncertainties arise from surface heat flux estimates and the representation of vertical mixing, but these are not expected to modify the relative importance or seasonal phasing of the dominant budget terms (Josey et al., 2013; Large et al., 1994).

Regarding vertical oxygen transport in the whole Levantine Basin, the weak upward transfer from the deep layer into the intermediate layer found in our results is consistent with the general scheme of circulation or oxygen cycle shown in previous studies (Mavropoulou et al., 2020; Powley et al., 2016; Roether and Schlitzer, 1991; Tanhua et al., 2013) describing a gradual upwelling of deep water originating from the Adriatic Sea or Aegean Sea. While a downward export of oxygen from the surface layer to the intermediate layer is simulated at the scale of the whole basin, the Rhodes Gyre exhibits an opposite pattern, with oxygen being transported upward from the intermediate layers to the surface. This upward input into the surface layer is balanced by a lateral export, particularly strong in winter, which takes place notably through the dispersal by subduction of the newly formed LIW at the periphery of the mixed patch, which is consistent with the observations reported by Malanotte-Rizzoli et al. (2003) for January 1995 during the POEM cruise and by Taillandier et al. (2022) during PERLE cruises.

Finally, our results on lateral oxygen exchanges are also in agreement with previous studies describing the general circulation in the Eastern Mediterranean Sea under BIOS cyclonic phases. During these cyclonic phases of the BIOS, the cyclonic Northern Ionian Gyre promotes an important eastward advection of Atlantic Water toward the Levantine Basin. As for the exchanges with the Ionian Sea, the general cyclonic circulation displays in the surface and intermediate layers an eastward inflow along the Libyan-Egyptian coast (Estournel et al., 2021). South of Crete, the flux reverses seasonally with an inflow from the Ionian in winter and an outflow in summer (Estournel et al., 2021). The net oxygen flux directed from the Ionian toward the Levantine Basin in surface and intermediate waters in this study results from a stronger southern input than the northern export at the Ionian-Levantine boundary. Regarding the exchanges with the Aegean Sea, a net outflow of LSW and LIW by the Asia Minor Current through the Cretan Straits was documented in several observational and modeling studies (Estournel et al., 2021; Millot and Taupier-Letage, 2005; Velaoras et al., 2014). During anticyclonic phases of the BIOS, the inflow from the Ionian Sea of Atlantic Water is reduced and significant changes in oxygen circulation may be expected and deserve further investigation.

In addition, extreme events may episodically modulate air–sea oxygen exchanges in the Levantine Basin. In particular, marine heatwaves, whose frequency and intensity have increased in the Eastern Mediterranean in recent decades (Aboelkhair et al., 2023; Darmaraki et al., 2024), can enhance upper-ocean stratification, reduce vertical ventilation, and decrease oxygen solubility through surface warming, potentially leading to transient oxygen anomalies in the upper and intermediate layers (Keeling et al., 2010; Schmidtko et al., 2017). Moreover, Medicanes, short-lived and spatially localized extreme events, could also impact biological dynamics and air–sea exchanges on their passage (Menna et al., 2023; Jangir et al., 2024, 2026; Reale et al., 2026). While a dedicated analysis of the impacts of those extreme events (marine heatwaves, Medicanes) is beyond the scope of the present study, their effects may contribute to short-term departures from the mean seasonal oxygen cycle and their integrated contribution to basin-scale and annual oxygen budget will need to be assessed in future works.

Beyond the Levantine Basin, the processes identified here have broader implications for regional biogeochemical dynamics and Earth system modelling. The quantified oxygen budget of LIW highlights the sensitivity of intermediate-water ventilation to circulation and atmospheric forcing at seasonal to interannual timescales in semi-enclosed basins (e.g. Tanhua et al., 2013; Schneider et al., 2014). Because LIW supplies the Eastern Intermediate Water and flows in the whole Mediterranean basin, variability in its oxygen content may propagate in the other water formation areas (South Aegean, South Adriatic, north-western Mediterranean), and as it constitutes the main precursor of the Mediterranean Overflow Water, can influence the oxygen and biogeochemical properties of intermediate waters beyond the Mediterranean in the Northeast Atlantic (Aldama-Campino and Döös, 2020; Stendardo and Gruber, 2012). Accurately representing these processes is therefore essential for regional biogeochemical models and Earth system models aiming to capture Mediterranean–Atlantic exchanges and their contribution to large-scale oxygen budgets.

5 Conclusion and future works

The study period was marked by contrasted atmospheric and hydrodynamic winter conditions. The confrontation of the model results with cruise and BGC-Argo float observations shows the capacity of the model to capture the general seasonal and spatial dissolved oxygen variability, as well as the main oxygen features in the Levantine Basin. These in situ observations, particularly from BGC-Argo floats and ship-based measurements, were essential for constraining and validating the simulations, without which the model outputs would not have reached their current level of reliability. The following conclusions can be drawn:

• The model results indicate a clear seasonal cycle for the oxygen air–sea flux. During winter, with the decrease in temperature, the increase in heat losses and intensified vertical mixing events, the surface layer is undersaturated in oxygen, resulting in atmospheric oxygen uptake. Undersaturation averaged over the whole basin reaches 2 % during winter. During the stratified period, primary production combined with the decrease of the temperature-dependent solubility in the thin mixed layer above the SOM leads to oxygen oversaturation and subsequent outgassing.

• At the annual scale, the Levantine Basin acts as a net sink for atmospheric oxygen, capturing 1.5 ± 0.6 mol O₂ m⁻² yr⁻¹ of oxygen. Most of the oxygen uptake occurs during winter, when it accounts for 10.7 ± 2.8 mol O₂ m⁻² yr⁻¹. The Rhodes Gyre absorbs atmospheric oxygen at a 2-fold higher rate than the entire Levantine Basin.

• Our budget shows that the surface layer of the Levantine Basin is a net source of dissolved oxygen for the intermediate waters, with winter vertical export of oxygen strongly modulated by the winter heat loss intensity. Regarding the exchanges with the surrounding seas, we found that oxygen is laterally transported into the Levantine Basin by surface and intermediate waters originating from the Ionian Sea. The lateral annual oxygen outflow toward the Aegean Sea is strongly enhanced by the heat loss intensity with exports 1.5 and 2.4 times higher during cold years in the surface and intermediate layer, respectively, compared to mild years.

• On an annual level, the Levantine Basin is found to act as a weak autotrophic ecosystem, with a net community production in the surface layer alternating between auto- and heterotrophic status influenced by the magnitude of the winter heat loss. In deeper depths and respiration resulted in an oxygen consumption of 2.0 ± 0.3 mol O₂ m⁻² yr⁻¹. Spatially, the Rhodes Gyre appears to be a major oxygen reservoir across the basin, contributing 41 % of the oxygen production of the whole surface layer.

This study represents a first step in our modeling of the dissolved oxygen dynamics in the Levantine Basin. While the quasi-permanent Rhodes Gyre dominates the basin-scale oxygen budget, transient cyclonic and anticyclonic mesoscale structures are also expected to contribute to oxygen variability outside this gyre on shorter time scales, in particular through nutrient upwelling or transport from nutrient-rich coastal waters offshore (Di Biagio et al., 2022; Pirro et al., 2024). Further investigations focusing on the specific role of these various cyclonic and anticyclonic eddies will be conducted in the future. Future work could also benefit from applying variability-based approaches, such as empirical orthogonal function (EOF) or regime-oriented analyses on model fields, to further disentangle the respective roles of physical and biogeochemical drivers across temporal scales (Di Biagio et al., 2023).

While the 7-year study period provides high-resolution insights into oxygen dynamics, it does not cover long-term climate shifts such as the EMT. Several studies suggest a decadal variability of dissolved oxygen across the whole water column linked to the dense water formations in the south Adriatic and Aegean seas and to the general eastern Mediterranean circulation, notably the reversal of the North Ionian Gyre (Ozer et al., 2020, 2022). Extended the simulation period, in addition to the implementation of a finer vertical resolution at key depths could contribute to examining this longer-term variability in the Levantine Basin and the connections between the sub-basins of the eastern Mediterranean. Improving the representation of intermediate water oxygen dynamics in the Mediterranean is also a necessary step toward better quantifying Mediterranean–Atlantic biogeochemical coupling and its sensitivity to future climate-driven changes in ventilation and circulation.