# JGR Oceans

## RESEARCH ARTICLE

10.1029/2022JC019245

**Key Points:**

* The seasonal atmospheric change is the main driver of destratification

* Winters with deep convection have below average levels of stratification that the atmospheric forcing has to overcome

* The Mistral winds have a low frequency signature that elevates the seasonal wind speeds in the winter

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

D. Keller Jr.,
dg.kllr.jr@gmail.com

**Citation:**

Keller, D., Jr., Givon, Y., Pennel, R., Raveh-Rubin, S., & Drobinski, P. (2024). Untangling the mistral and seasonal atmospheric forcing driving deep convection in the Gulf of Lion: 1993–2013. *Journal of Geophysical Research: Oceans*, 129, e2022JC019245. https://doi.org/10.1029/2022JC019245

Received 27 AUG 2022
Accepted 29 MAY 2024

**Author Contributions:**

**Conceptualization:** Shira Raveh-Rubin, Philippe Drobinski

**Formal analysis:** Douglas Keller Jr., Yonatan Givon

**Investigation:** Philippe Drobinski

**Methodology:** Douglas Keller Jr., Yonatan Givon, Romain Pennel

**Project administration:** Shira Raveh-Rubin, Philippe Drobinski

**Supervision:** Romain Pennel, Shira Raveh-Rubin, Philippe Drobinski

**Writing – original draft:** Douglas Keller Jr.

**Writing – review & editing:** Douglas Keller Jr., Yonatan Givon, Romain Pennel, Shira Raveh-Rubin, Philippe Drobinski

© 2024 The Author(s).
This is an open access article under the terms of the Creative Commons Attribution-NonCommercial License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited and is not used for commercial purposes.

<figure><img src="image_1.png" /></figure>

# Untangling the Mistral and Seasonal Atmospheric Forcing Driving Deep Convection in the Gulf of Lion: 1993–2013

Douglas Keller Jr.¹, Yonatan Givon², Romain Pennel¹, Shira Raveh-Rubin², and Philippe Drobinski¹

¹LMD/IPSL, École Polytechnique, Institut Polytechnique de Paris, ENS, PSL Research University, Sorbonne Université, CNRS, Palaiseau, France, ²Department of Earth and Planetary Sciences, Weizmann Institute of Science, Rehovot, Israel

**Abstract** Deep convection occurs periodically in the Gulf of Lion, in the northwestern Mediterranean Sea, driven by the seasonal atmospheric change and Mistral winds. To determine the variability and drivers of both forcings, multiple 1 year ocean simulations were run, spanning from 1993 to 2013. Two sets of simulations were performed: a control and seasonal set, the first forced by unfiltered atmospheric forcing and the other by filtered forcing. The filtered forcing was bandpass filtered, retaining the seasonal and intraday aspects but removing the high frequency phenomena. Comparing the two sets allows for distinguishing the effects of the high frequency component of the Mistral on the ocean response. During the preconditioning phase, the seasonal forcing was found to be the main destratifying process, removing on average 46% of the stratification needed for deep convection to occur, versus the 28% removed by the Mistral. Despite this, each forcing triggered deep convection in roughly half of the deep-convection events. Sensible and latent heat fluxes were found to be the main drivers of the seasonal forcing during deep-convection years, removing 0.17 and 0.43 m²s⁻² of stratification, respectively. They were themselves driven by increased wind speeds, believed to be the low frequency signal of the Mistral, as more Mistral events occur during deep-convection winters (34% vs. 29% of the preconditioning period days). An evolving seasonal forcing in a changing climate may have significant effects on the future deep convection cycle of the western Mediterranean Sea.

**Plain Language Summary** Deep convection occurs periodically in the Gulf of Lion (located in the northwestern Mediterranean Sea), when water at the surface of the ocean is cooled enough to mix freely with the deeper water below, sometimes reaching the sea floor. It's an important part of the overall circulation of the Mediterranean Sea that leads to an explosion in the phytoplankton population in the following spring. In the gulf, the surface cooling is caused by the atmospheric transition from summer to winter and the Mistral winds. The latter is a cool, dry northerly wind that flows through the Rhône Valley out over the gulf. In our study, we ran ocean simulations that included and excluded the non-seasonal effects of the Mistral to determine the importance of the seasonal and Mistral forcing on deep convection. The seasonal forcing was found to have a larger role, and contained a low frequency part of the Mistral, elevating the average wind speeds during the winter. Changes in the seasonal forcing and ocean water composition will need to be studied to understand the evolution of deep convection in the Gulf of Lion and its consequences on the Mediterranean Sea dynamics and biology in a changing climate.

## 1. Introduction

Deep convection, or open-ocean convection, occurs in the higher latitude regions of the world and is an important ocean circulation process (Marshall & Schott, 1999). It is formed when the stable density gradient along the ocean column is eroded by surface buoyancy loss, leading to an overturning that can span the entire depth of the water column. In the western basin of the Mediterranean Sea, this process can occur in the Gulf of Lion and assists in the thermohaline circulation of the sea (Robinson et al., 2001) by forming the Western Mediterranean Deep Water (WMDW). When it does occur, the WMDW produced spreads out along the lower layers and bottom of the northwest basin (MEDOC, 1970). Figure 13 in Testor et al. (2018) illustrates this process, additionally showing the submesoscale coherent vorticies that are produced and aid in the lateral transport (Bosse et al., 2016). Some of the WMDW is transported along the northern boundary current toward the Balearic Islands (Send & Testor, 2017), and some of it completes the general circulation by flowing down toward the Algerian Basin and the Strait of Gibraltar (Beuvier et al., 2012; Testor & Gascard, 2003). In the Gulf of Lion, deep convection also plays an important role in the marine biology of the region, as the springs following deep convection events also
experience increased phytoplankton blooming (Severin et al., 2017), due to the increased levels of nutrients and oxygenation from the mixing process (Coppola et al., 2017; Severin et al., 2017).

Significant deep convection events occur every few years in the Gulf of Lion (Bosse et al., 2021; Houpert et al., 2016; Marshall & Schott, 1999; Mertens & Schott, 1998; Somot et al., 2016), driven by the Mistral and Tramontane winds. These sister, northerly flows bring cool, continental air through the Rhône Valley (Mistral) and the Aude Valley (Tramontane), leading to large heat transfer events with the warmer ocean surface (Drobinski et al., 2017; Flamant, 2003). These large cooling, evaporative events destabilize the water column in the Gulf of Lion, and are a primary source of buoyancy loss leading to deep convection (Herrmann et al., 2010; Lebeaupin-Brossier & Drobinski, 2009; Madec, Chartier, & Crépon, 1991; Madec, Delecluse, et al., 1991). The other main source of buoyancy loss in the region is the seasonal atmospheric change and reduction of solar heating during the winter (Keller et al., 2022).

Multiple studies have looked at this region and phenomenon. Some recent examples include Waldman et al. (2018) who looked at the inherent chaotic variability of deep convection in the Gulf of Lion and found that the process is random 18% of the years from 1979 to 2013, with the most intrinsic variability in the convection patch, away from the coast. Herrmann et al. (2017) proposed and analyzed a method using satellite altimetry and chlorophyll data to estimate dense water formation in the gulf over roughly two decades. Another work, Amitai et al. (2021), investigated the connection between the source of intermediate and deep water formation in the eastern Mediterranean on the Gulf of Lion deep convection. They found that the Eastern Mediterranean Transient (EMT) could be modeled as a shift between two dynamically stable regimes in the eastern part of the sea, one regime featuring the primary source of dense water formation in the east in the Adriatic Sea and the other featuring it in the Aegean Sea, which can either impede or enhance deep convection in the Gulf of Lion. Past works have also investigated the effects of spatial resolution and temporal frequency of atmospheric forcing on the phenomenon in the gulf, such as Herrmann and Somot (2008) and Madec, Chartier, and Crépon (1991), respectively. However, separating the effects of atmospheric forcing by frequency over multiple years has yet to be investigated.

The annual stratification cycle of the Gulf of Lion regulates the occurrence of deep convection events. It comprises of a destratification phase and restratification phase that is roughly sinusoidal in appearance. These two phases form due to the net heat flux into the ocean surface changing sign roughly at the spring and fall equinoxes: positive between March and September and negative between September and March. When the net heat flux is positive, the ocean column is being heated, increasing its stability, hence an increase in stratification from March to September. When the net heat flux is negative, the ocean column is being cooled, reducing its stability, thereby decreasing its stratification from September and March. The net heat flux gains its shape from its four main components: solar heating, infrared cooling, the sensible heat flux, and the latent heat flux. The solar heating gives the net heat flux its sinusoidal shape. The infrared cooling, sensible heat flux, and latent heat flux shift this sinusoidal shape negative, causing it to flip sign at the spring and fall equinoxes. The asymmetries in the net heat flux come from the sensible and latent heat fluxes, causing the sinusoidal shape to be distorted slightly in the winter (Keller et al., 2022).

The aforementioned destratification phase is also known as the preconditioning phase, as the fall and winter surface heat fluxes precondition the ocean vertical column for deep convection in the spring (Keller et al., 2022; Marshall & Schott, 1999; MEDOC, 1970). However, in addition to these surface fluxes, lateral advection within the ocean column can have an impact on the vertical stability. For example, Margirier et al. (2020) showed the effects of a warming and more saline Levantine Intermediate Water (LIW) in the gulf. The warmer and more saline LIW increased the vertical stability, causing the ocean vertical column to require more cooling fluxes during the fall, winter, and early spring to overcome the increase in stratification.

If the cooling from the sensible and latent heat fluxes is large enough during the preconditioning phase (the infrared cooling tends to remain constant as it depends on the sea surface temperature), then a third phase appears: the deep convection phase. This occurs when the sensible and latent heat fluxes reduce the stratification to point it can overturn. These three phases then form the canonical deep convection cycle (MEDOC, 1970; The Lab Sea Group, 1998). For this study, we are focusing on the destratification (preconditioning) and deep convection phases, as they drive the variability of this cycle in the Gulf of Lion.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** The domains of both the WRF domain from the RegIPSL coupled WRF/ORCHIDEE simulation within the Med-CORDEX framework, (a), and the NEMOMED12 configuration domain, (b). The region of interest, the NW Med., is outlined by the box. This region is later used in Figure 4.</figcaption></figure>

In a prior, closely related paper, Keller et al. (2022), we determined the importance of the seasonal atmospheric change with regards to its impact on the destratification phase and discovered it was a more significant source of destratification than the Mistral/Tramontane wind events (referred to as just Mistral events), providing roughly 2/3 of the destratification for the 2012–2013 winter. The current study continues this investigation and looks into the variability of the contribution to destratification for each component, the seasonal and the Mistral, over multiple years. 20 years of the Mediterranean Sea, from 1 July 1993 to 30 June 2013, were simulated using the NEMO ocean model. NEMO was driven by two sets of WRF/ORCHIDEE atmospheric data: a control set and a filtered (seasonal) set. This resulted in two sets of the simulated ocean data: one set including the effects of the Mistral and the seasonal effects, and the other set just including the seasonal effects, allowing us to separate the effects due to the high frequency component of the Mistral.

In particular, our findings determine the variability of both the seasonal and Mistral based contributions to destratification, and the primary components and their drivers, of the seasonal contribution leading to deep convection. Our findings also address questions posed by Keller et al. (2022) that were outside the scope of that study. These questions can be summarized as the following:

1. Does the Mistral trigger deep convection, or does the seasonal atmospheric change trigger it?

2. Does the maximum seasonal stratification play a role in deep convection?

3. Does the timing of the seasonal stratification minimum matter and can the Mistral contribution overcome a restratifying seasonal stratification?

The paper is organized in the following way. The NEMO model, atmospheric forcing data and filtering, additional methodology, and observational data are described in the Methodology section (Section 2). The ocean model outputs are validated using the observational data in the Model Validation section (Section 3). The results of the seasonal and Mistral contributions are presented and discussed in the Results and Discussion section (Section 4), along with addressing the questions posed by Keller et al. (2022), followed by the conclusion (Section 5).

## 2. Methodology

To separate the effect of the high frequency component of the Mistral and seasonal aspects of the atmospheric forcing, two sets of ocean simulations simulating the Mediterranean Sea were carried out: one control set and one seasonal set. The seasonal set had part of its atmospheric forcing filtered to remove the Mistral from the forcing, thereby allowing the differences between the two ocean simulation sets to reflect the effect the Mistral has on the ocean (the filtering process is described in Section 2.2). As the Mistral is the main intra-monthly phenomenon that occurs during the winter in the Gulf of Lion (Givon et al., 2021; Keller et al., 2022), the seasonal ocean simulations reflect the ocean response just due to the seasonal atmospheric changes in the region. The two ocean simulation sets are performed on a per year basis from the same initial conditions pulled from the MEDRYS reanalysis (Hamon et al., 2016). For example, one control and seasonal simulation pair was run from 1 July 1993 to 30 June 1994, both initialized by the data from 1 July 1993 from MEDRYS. The same was performed from 1 July 1994 to 30 June 1995, again reinitialized on July 1st from the data from 1 July 1994 from MEDRYS, and so on, until 30 June 2013. This was done to allow for the assumption that processes outside the northwestern Mediterranean subdomain in Figure 1b that are affected by the filtering, such as the formation of TIW (Tyrrhenian Intermediate Water), which is driven by winter surface heat fluxes (Napolitano et al., 2019), have a negligible
impact on the Gulf of Lion processes during the comparison of per year ocean simulations. This assumption is corroborated by the slow movement of intermediate and dense water, which is on the order of 2 years for intermediate waters to travel from the Strait of Sicily to the Gulf of Lion (Amitai et al., 2021), through the Tyrrenian Sea. It is also corroborated by the roughly year time scale for newly formed WMDW (Western Mediterranean Deep Water) to move into the southern Algerian Basin (Beuvier et al., 2012) and the order of decades time scale for total circulation of the Mediterranean Sea (Millot & Taupier-Letage, 2005).

## 2.1. NEMO

The Nucleus for European Modeling of the Ocean (NEMO) ocean model (https://www.nemo-ocean.eu/; last accessed: 16 August 2022; (Madec, 2008; Madec et al., 2022)) was used to simulate the Mediterranean Sea in 1 year runs for 20 years, as described above, from 1 July 1993 to 30 June 2013. The model was run in bulk configuration, utilizing the following parameterized equations:

$$Q_E = \rho_{a,0} \Lambda C_E (\Delta q) |\Delta \vec{u}| \quad (1)$$

$$Q_H = \rho_{a,0} c_p C_H (\Delta \theta) |\Delta \vec{u}| \quad (2)$$

$$Q_{LW} = Q_{LW,a} - \epsilon \sigma SST_K^4 \quad (3)$$

$$\tau = \rho_{a,0} C_D \Delta \vec{u} |\Delta \vec{u}| \quad (4)$$

Where $Q_E$, $Q_H$, $Q_{LW}$, and $\tau$ are the latent heat, sensible heat, longwave radiation fluxes and the surface shear stress, respectively. $z$ is the height above the sea surface where the atmospheric variables are provided at (at 2 m for temperature and humidity and at 10 m for the wind velocity), with the naught values ($_0$) at the sea surface. $\vec{u}$ is the horizontal wind vector, with $\Delta \vec{u} = \vec{u}_z - \vec{u}_0$ as the difference between the wind velocity and sea surface current. $\Delta q = q_z - q_0$ and $\Delta \theta = \theta_z - \theta_0$; $q$ and $\theta$ are the specific humidity and potential temperature of air, respectively. $\Lambda$ and $c_p$ are the latent heat of evaporation and the specific heat of water, respectively. $\rho_a$ is the density of air. $SST_K$ is the sea surface absolute temperature. $\epsilon$ is the sea surface emissivity, $\sigma$ is the Stefan-Boltzmann constant, and $Q_{LW,a}$ is the atmospheric longwave radiation. The coefficients $C_E$, $C_H$, and $C_D$ are the parameterized coefficients of latent heat, sensible heat, and drag, respectively, and are defined in Large and Yeager (2004, 2008).

$Q_{net}$, the net downward heat flux, is the summation of the components in the following equation ((Estournel et al., 2016; Large & Yeager, 2004); ignoring snowfall):

$$Q_{net} = Q_{SW} + Q_{LW} + Q_H + Q_E \quad (5)$$

Where $Q_{SW}$ is the downward shortwave radiation.

The NEMO model was also run in the NEMOMED12 configuration, using NEMO v3.6. NEMOMED12 is described, with boundary conditions, in Lebeaupin-Brossier et al. (2011), Beuvier et al. (2012), Hamon et al. (2016), and Waldman, Somot, et al. (2017). This configuration has been used extensively to study dense water formation around the Mediterranean, such as by Waldman, Somot, et al. (2017), Waldman, Herrmann, et al. (2017), and Keller et al. (2022) to investigate the 2012–2013 deep convection event and by Waldman et al. (2018) to study the chaotic variability of deep convection in the Gulf of Lion. A similar, coarser configuration, called NEMOMED8 was used by Herrmann et al. (2010) and Somot et al. (2016) to study the inter-annual variability of deep convection and by Adloff et al. (2015) to see the predicted changes in the Mediterranean due to climate change.

A brief description of the configuration NEMOMED12 follows: the domain covers the Mediterranean Sea and a portion of the Atlantic Ocean (see Figure 1b). The horizontal resolution is at 1/12°, with a horizontal distance resolution of about 7 km, coming from the ORCA 1/12° global ocean grid (Hurlburt et al., 2009). There are 75 vertical levels ranging from 0.5 m thick at the surface to 235 m thick at the bottom. The simulations were run at 720 s time steps. The Atlantic buffer zone is used to represent the exchanges between the two bodies of water at the Strait of Gibraltar, and its sea surface height (SSH) fields are restored toward the ORAS4 global ocean reanalysis (Balmaseda et al., 2013), with a no-slip lateral boundary condition applied. The 3-D temperature and
salinity fields of the buffer zone are restored toward the MEDRYS reanalysis (Hamon et al., 2016). The Black Sea, runoff of 33 major rivers, and coastal runoff are represented by climatological data from Ludwig et al. (2009) and are provided as monthly values. The initial conditions for each 1 year run were pulled from the MEDRYS reanalysis (Hamon et al., 2016) as mentioned before.

## 2.2. Atmospheric Forcing

The atmospheric forcing used in the simulations were the output of a RegIPSL simulation, the regional climate model of IPSL (Guion et al., 2021b) (https://gitlab.in2p3.fr/ipsl/lmd/intro/regipsl/regipsl; last accessed: 26 August 2022), which used the coupling of the Weather Research and Forecasting Model (WRF) (Skamarock et al., 2008) and the ORCHIDEE Land Surface Model (Krinner et al., 2005). The run is a hindcast simulation (ERA Interim downscaling), performed at 20 km resolution, spanning the period of 1979–2016, within the HyMeX (Drobinski et al., 2014) and Med-CORDEX framework (Ruti et al., 2016). The u and v wind components, specific humidity, potential temperature, shortwave and longwave downward radiation, precipitation, and snowfall at 3 hourly temporal resolution were all used to force the NEMO ocean simulations.

For the control simulation set, the forcing were used as is. For the seasonal simulation set, the u and v wind components, specific humidity, and potential temperature were filtered over the entire domain of the WRF forcing (Figure 1a). These variables are the primary variables that have the imprint of the Mistral in them and are involved with surface flux calculations utilizing the bulk formulae (Equation 4). Thus filtering them is equivalent to filtering out the Mistral. The variables relating to the radiation and precipitation fluxes were left unchanged, as they are not significantly altered by the Mistral (and therefore do not contain its signal). The filtering process was performed by a moving window average:

$$ \chi_i = \frac{1}{2N + 1} \sum_{j=i-N}^{i+N} x_j \qquad (6) $$

Where $\chi_i$ is the averaged (filtered) value at index $i$ of a time series of variable $x$ with length $n$, where $i = 0 \rightarrow n$. The window size is equal to $2N + 1$, which, in this case, is equal to 31 days. The ends have a reduced window size for averaging, and thus show edge effects. The edge effects did not affect the forcing used for the NEMO simulations, however, as they were before and after the overall ocean simulation beginning and end dates.

The moving window average was applied to each time point per day over a 31 days window. that is, for 3 hourly data, the time series is split into 8 separate series, one for each timestamp per day, (00:00, 03:00, 06:00, etc.) and then each series is averaged with a moving window. The eight window averaged series are then recombined into a single time series. This was done to retain the intra-day variability, yet smooth the intra-monthly patterns, as the diurnal cycle has been shown to retard destratification by temporarily reforming a stratified layer at the sea surface during slight daytime warming. This diurnal restratification has to be overcome first before additional destratification of the water column can continue during the next day (Lebeaupin-Brossier et al., 2011, 2012) and is shorter than a typical Mistral event length of a little over 5 days (Keller et al., 2022). An example of the filtering can be seen in Figure 2 of Keller Jr. et al. (2022). The filtering removes the short term, anomaly scale forcing from the forcing data set (the phenomena with under a month timescale), effectively removing the Mistral's influence on the ocean response. This creates two separate forcing data sets: one with the anomaly forcing included, attributed to the Mistral and hence called the Mistral forcing, and one with just the seasonal forcing, leading to the designation of “control” and “seasonal” for the unfiltered and filtered data sets, respectively. The resulting control and seasonal forcing fluxes are shown in Figures S1 and S2 in Supporting Information S1.

The main assumption of performing this filtering is the Mistral primarily acts on the short term, anomaly scale forcing (high frequency forcing). This found to be a fairly effective assumption when separating Mistral and seasonal effects in Keller et al. (2022). However, there is a seasonal component to the Mistral forcing that is not removed with this filtering. Mistral events occur more frequently in the winter than in summer (Givon et al., 2021), which appears in lower frequencies of the atmospheric forcing. This will be discussed more in Section 4.1.1.
## 2.3. Stratification Index

The stratification index, SI, is a useful measure of the stability of the ocean column (equivalently known as “convection resistance” (Frajka-Williams et al., 2014)). It builds from the non-penetrative growth of the mixed layer, a reasonable assumption for the ocean mixed layer (Keller et al., 2022; Somot, 2005; Turner, 1973). It compresses the Brunt-Väisälä frequency, $N^2$, over the depth of the water column into a single index:

$$ SI = \int_{0}^{D} N^{2} zdz \qquad (7) $$

Where $z$ is the depth and $D$ is the depth of ocean column. As $N^2$ is proportional to the vertical density gradient, SI provides a 0 dimensional metric to measure the stratification of the ocean column. We will use it to track the stratification of the Gulf of Lion and the occurrences of deep convection. Consequently, the stratification indexes from the control and seasonal simulation sets are SI and $SI_S$, respectively, with the difference, $\delta SI = SI - SI_S$, being the stratification induced by the Mistral.

## 2.4. Simple Model

To separate the different seasonal drivers of deep convection, a simple model that relates the seasonal stratification index, $SI_S$, to the seasonal net surface heat flux, $Q_{net,S}$ is used (Keller et al. (2022) Equation 16):

$$ \frac{dSI_S}{dt} = \frac{g}{2\rho c_p T_0} Q_{net,S} \approx 10^{-9} \times Q_{net,S} \qquad (8) $$

Where $c_p$ is the specific heat capacity of water (taken as 4,184 Jkg$^{-1}$K$^{-1}$), $g$ is gravity (taken as 9.81 $m^{-1}$), $\rho$ is the reference density of water (taken as 1,000 kgm$^{-3}$), and $T_0$ is the reference temperature (taken as 290 K; the average seasonal sea surface temperature over the 20 year period). Utilizing these values, $\frac{g}{2\rho c_p T_0} \approx 10^{-9} m^4 J s^{-2}$.

$Q_{net,S}$ can be further separated into its individual, $i$, components through Equation 5, allowing us to estimate the components' individual contribution to destratification/restratification by integrating over a selected interval of time ($t_0-t_1$):

$$ SI_i = 10^{-9} \times \int_{t_0}^{t_1} Q_i dt \qquad (9) $$

Note, however, that Equation 9 neglects the advection terms on the low frequency time scale, and assumes that the stratification in the gulf is driven completely by surface fluxes. This isn't entirely accurate, as Herrmann et al. (2010) has shown non-insignificant lateral fluxes in the salt and heat content exist, however at roughly a quarter of the magnitude of the surface fluxes.

Similarly, the effect of a Mistral event, $k$, on destratification, $\delta SI_k$, can be calculated using Equation 17 from Keller Jr. et al. (2022):

$$ \Delta\delta SI_k = \delta SI_k(t_k + \Delta t_k) - \delta SI_{k-1}(t_k) = \left[ \delta SI_{k-1}(t_k) + \frac{D^2}{2} \frac{\delta F_k}{\alpha_d} \right] (e^{-\alpha_d \Delta t_k} - 1) \qquad (10) $$

This solution assumes the Mistral behaves as a pulse function, shown in Figure 2, with the width of the pulse as the Mistral event duration, $\Delta t_k$, its height as the Mistral event's “strength” or its cooling effect from the surface heat flux driven by the Mistral, $\delta F_k$, and the distance between pulse starts, $t_k$, as the period between of multiple events, $\Delta \tau_k$. $\alpha_d$ is the restoration coefficient, which is effectively the horizontal gradient of $\delta SI$ during the event (Keller et al., 2022), which comes about from the effects of advection. We set $\alpha_d$ equal to 0.23 day$^{-1}$, which is taken from Keller et al. (2022). The calculation for $\delta F_k$ is provided in Appendix A with the resulting values in Tables S1 in Supporting Information S1 through Table S9 in Supporting Information S1.

As an aside, Equations 8 and 10 are both built from the energy equation for incompressible flow (derivation in Keller et al. (2022)). Within that derivation, it is assumed that Brunt-Väisälä frequency is constant throughout the water column, thereby assuming the density gradient is constant in the vertical. In reality, the density gradient is
<figure><img src="image_3.png" /><figcaption>**Figure 2.** An example time series of pulse functions used to analytically represent Mistral events. Each event starts at time $t_k$, has a duration of length $\Delta t_k$, and period before the beginning of the following event, $\Delta \tau_k$. $\delta F_k$ represents the event's "strength" or cooling effect due the surface heat fluxes driven by the Mistral.</figcaption></figure>

not constant, and experiences sharp increases, particularly at the boundaries of different water masses. Despite this, the equations are still useful to explain and breakdown the different contributors to changes in stratification, albeit with this caveat in mind.

## 2.5. Mistral Events

To use Equation 10 and investigate more into how the Mistral affects the Gulf of Lion, the Mistral events during the 1 July 1993–30 June 2013 period are pulled from the Mistral data set provided by Keller Jr. et al. (2022) and Givon et al. (2021) (https://medcyclones.eu/data/; last accessed 23 August 2022). Mistral events are essentially determined by the simultaneous presence of northerly flows in the Rhône Valley and over the Gulf of Lion, in conjunction with a low pressure system over the Ligurian Sea. More details are provided in Keller et al. (2022) and Givon et al. (2021).

## 2.6. Argo and CTD Profiles

To validate the control set of the ocean simulations, Argo and CTD vertical profile observations from the period of 1 July 1993 to 30 June 2013 were retrieved from the Coriolis database (https://www.coriolis.eu.org/Data-Products/Data-selection; last accessed: 23 Aug. 2022; some of which were provided by the Atlantic Oceanographic and Meteorological Laboratory). 2929 temperature and salinity in situ profiles were retrieved. 1949 profiles were from Argo profiling and 980 were from CTD profiling (vertical distribution breakdown in Table 1 and locations in Figure S3 in Supporting Information S1). These vertical profiles were compared to the model vertical profiles to determine and verify the accuracy of the model. The model outputs salinity in terms of practical salinity, in units of PSU, which is the same as the observational data. However, for temperature, the model outputs potential

<caption>**Table 1**</caption>

Number of and Start and end Dates for the Argo and CTD Profiles From the Coriolis Database for the Period 1 July 1993 to 30 June 2013

<table><thead><tr><th></th><th>Start</th><th>End</th><th>0–150 m</th><th>150–600 m</th><th>Below 600 m</th><th>Total</th></tr></thead><tbody><tr><td>Argo</td><td>2005-01-01</td><td>2011-12-30</td><td>1945</td><td>1684</td><td>1479</td><td>1949</td></tr><tr><td>CTD</td><td>1993-07-05</td><td>2013-06-29</td><td>978</td><td>481</td><td>205</td><td>980</td></tr><tr><td>Total</td><td>—</td><td>—</td><td>2923</td><td>2165</td><td>1684</td><td>2929</td></tr></tbody></table>

*Note.* The number of profiles used for the spatial distribution of bias in the layers from the surface to 150 m, from 150 to 600 m and below 600 m in depth are shown in their respective columns. These layers roughly correspond the typical depths of the Atlantic Water, Levantine Intermediate Water, and West Mediterranean Deep Water masses, respectively (Houpert et al., 2016; Margirier et al., 2020).
temperature, whereas the observed temperature is provided in terms of in situ temperature measurements. To make a direct comparison, the observational temperature data was converted to potential temperature with the GSW-Python python package (Firing et al., 2021), which uses the TEOS-10 ocean equation of state for the conversion (https://www.teos-10.org/index.htm; last accessed 23 August 2022).

## 2.7. Argo Mixed Layer Depths

In addition to the Argo and CTD profiles, the control set of the ocean simulation results for the mixed layer depth (MLD) were also compared to MLDs calculated from the global Argo data set (separate from our set in the above section; the MLD calculations were performed by Holte et al. (2017a) and provided at https://mixedlayer.ucsd.edu; last accessed: 22 February 2023) (Holte et al., 2017a). From this data set, 3846 MLD calculation points were used for our time period and domain (the area covered by Figures 4–6). For this data set, the MLD was determined by a density based algorithm (Holte & Talley, 2009) (also see http://mixedlayer.ucsd.edu for more details; last accessed: 22 February 2023). Their algorithm is fairly in-depth and is only briefly described here via a quote from Holte and Talley (2009): “the algorithm models the profile’s general shape, calculates a suite of possible MLD values, and then looks for groupings and patterns within the possible MLDs to select the final MLD estimate for each profile.”

## 2.8. QuikSCAT Winds

To validate the control set of the wind forcing, the QuikSCAT gridded surface wind data set was utilized (https://cersat.ifremer.fr/Data/Latest-products/Reprocessing-Daily-QuikSCAT-Surface-Wind-Fields; last accessed: 22 February 2023). The data set provides daily surface wind speeds (meridional and zonal) collected by the QuikSCAT satellite scatterometer. Wind speeds are available from July 1999 to November 2009 at 0.25° spatial resolution and will be compared over its entire duration to the unmodified (control) 10 m wind speeds used to force the ocean model.

# 3. Model Validation

## 3.1. Temperature and Salinity

Each profile of calculated potential temperature (see Section 2.6) and salinity was compared to the model profile from the nearest grid point in the NEMOMED12 grid and nearest time stamp (daily temporal resolution for the ocean simulation data; the model data was interpolated vertically to match the levels of the observations). The bias (model minus observation) and root mean squared error (RMSE) were calculated from the comparisons.

To look at the vertical distribution of bias and RMSE, the observations and nearest model data were vertically binned (55 bins) according to depth. The bias was then calculated per observation/model result pair. The mean and standard deviation of the bias per each bin are plotted in Figures 3a and 3b, for potential temperature and salinity, respectively. For each bin, the RMSE was computed, and is shown in Figures 3c and 3d, for potential temperature and salinity, respectively.

We additionally calculated the resulting density bias due to these biases via the following equation (available in the NEMO model (Madec et al., 2022)):

$$ \rho = \underbrace{-a_0(1 + 0.5\lambda_1 T_a + \mu_1 z)T_a}_{\rho_T} + \underbrace{b_0(1 - 0.5\lambda_2 S_a - \mu_2 z)S_a}_{\rho_S} - \nu T_a S_a + \rho_0 \quad (11) $$

$\rho_T$ and $\rho_S$ are the strictly temperature and salinity components of the equation. The remaining portion of the equation is not investigated in terms of bias. It is relatively small compared to the prior components and is therefore neglected. $T_a = T - 10 (°C)$ and $S_a = S - 35 (PSU)$, where T and S are the temperature and salinity of the ocean. In this equation, there are no differences between conservative and potential temperature nor practical and absolute salinity. $\rho_0$ is equal to 1,026 kg/m³. Table 2 holds the coefficients not directly named. The resulting calculated biases in $\rho$ are shown in plots (e) and (f) of Figure 3.

As seen in Figure 3, most of the differences between the model and observations lie within the first few hundred meters of the ocean column. The largest differences and variability in the bias (between the control simulation set
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Vertical distribution of bias (model minus observation) and RMSE from the comparison of our control set model results and combined Argo/CTD observations. (a) and (b) show the mean and standard deviation of the bias for potential temperature and salinity, respectively. The mean is the solid red line, with the shading representing the are encompassed by ±1 standard deviation. (c) and (d) show the RMSE for potential temperature and salinity, respectively. (e) and (f) show the density bias induced by the biases in (a) and (b), for temperature and salinity, respectively, calculated with Equation 11.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Spatial distribution of the bias from the comparison of our control set model results and combined Argo/CTD observations. The top row shows the bias in the potential temperature for the three water layers: 0–150, 150–600, and below 600 m for (a)–(c), respectively. The bottom row shows the same for salinity. The black box from 42 to 42.5°N and from 4.25 to 5°E bounds the spatial averaging performed in Section 4. For visualization purposes, the bias has been linearly interpolated via triangulation to the NEMO model grid points.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Spatial bias of the mixed layer depth (MLD) between the model and Argo observations, as well as the model MLD and *SI*. (a) shows the difference between the model and Argo MLD data, in meters (the spatial visualization is linearly interpolated like in Figure 4), and (b) shows the average yearly maximum MLD (in colored contours), with the annual *SI* minima overlaid as “x”s. (c) shows the corresponding temporal monthly mean and variability (in standard deviation) of the MLD bias. In this figure, the MLD is taken a positive number, therefore a positive bias means the model overestimated the MLD relative to the Argo observations. The black bounding box first mentioned in Section 3.1 and first shown in Figure 4 is provided for reference.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Spatial distribution of the bias between the 10 m magnitude of wind speed of the control set of the WRF-ORCHIDEE wind forcing and the QuikSCAT surface wind speed magnitudes. The subplots (a)–(c) show the bias, RMSE, and the difference in 95th quantile (model minus observations), respectively. (d) shows the correlation of the model and observations. Again, the black bounding box first mentioned in Section 3.1 is provided for reference.</figcaption></figure>

results and the observations) are found at the surface, with a mean bias and RMSE of +0.40°C and 1.18°C, for potential temperature, and -0.04 PSU and 0.01 PSU, for salinity. Below 600 m, the bias and RMSE are much smaller, with the mean bias and RMSE averaging at -0.005°C and 0.08°C, for potential temperature, and +0.004 PSU and 0.021 PSU, for salinity. This is reflected in the resulting larger, negative density biases of -0.083 kg/m³ due to temperature and -0.033 kg/m³ due to salinity at the surface and smaller, by an order of magnitude, density biases of 0.002 kg/m³ due to temperature and 0.003 kg/m³ due to salinity below 600 m. Most of the variability is in the upper layers.

The larger variabilities and differences in the upper layers can be explained by the high variability in the representation of the surface fluxes. For example, the diurnal cycle isn't captured in the daily temporal resolution of the ocean model data. The sea surface layer destratifies and restratifies with the diurnal cycle, as noted by Lebeaupin-Brossier et al. (2011, 2012), whereas the lower layers are less effected, hence showing less error

between the observations and model output. Additionally, the sea surface fluxes are turbulent in nature, and aren't always exactly represented, again leading to biases in the uppermost layers of the ocean. These differences in the representation can lead to changes in the mixed layer depth, particularly during periods of large vertical mixing, and differences in the representation of the mesoscale eddies formed during and after vertical mixing events, which can lead to spatial and temporal discrepancies between the model and observations due to their turbulent nature. With these caveats noted, the model is fairly representative of the vertical column in the Gulf of Lion, with slightly warmer and fresher surface waters relative to observations and fairly accurate temperature and salinity for the deeper waters.

To see if there is any notable features in the spatial distribution of bias, the averaged biases of the water masses in three different layers were examined: 0 to 150, 150 to 600, and below 600 m. These layers roughly correspond the

<table><caption>**Table 2**

Coefficients for Equation 11, Available in the NEMO Ocean Model (Madec et al., 2022)</caption><thead><tr><th>Coeff.</th><th>Value</th><th>Units</th><th>Desc.</th></tr></thead><tbody><tr><td>a<sub>0</sub></td><td>1.6550 × 10<sup>-1</sup></td><td>kgm<sup>-3</sup>°C<sup>-1</sup></td><td>Linear thermal exp. coeff.</td></tr><tr><td>b<sub>0</sub></td><td>7.6554 × 10<sup>-1</sup></td><td>kgm<sup>-3</sup>PSU<sup>-1</sup></td><td>Linear haline exp. coeff.</td></tr><tr><td>λ<sub>1</sub></td><td>5.9520 × 10<sup>-2</sup></td><td>°C<sup>-1</sup></td><td>Cabbeling coeff. for T<sup>2</sup></td></tr><tr><td>λ<sub>2</sub></td><td>5.4914 × 10<sup>-4</sup></td><td>PSU<sup>-1</sup></td><td>Cabbeling coeff. for S<sup>2</sup></td></tr><tr><td>ν</td><td>2.4341 × 10<sup>-3</sup></td><td>°C<sup>-1</sup>PSU<sup>-1</sup></td><td>Cabbeling coeff. for TS</td></tr><tr><td>μ<sub>1</sub></td><td>1.4970 × 10<sup>-4</sup></td><td>m<sup>-1</sup></td><td>Thermobaric coeff. for T</td></tr><tr><td>μ<sub>2</sub></td><td>1.1090 × 10<sup>-5</sup></td><td>m<sup>-1</sup></td><td>Thermobaric coeff. for S</td></tr></tbody></table>
typical depths of the surface waters, the LIW, and WMDW, respectively (Houpert et al., 2016; Margirier et al., 2020). They are plotted in Figure 4, with subplots (a), (b), and (c) showing temperature and (d), (e), and (f) showing salinity, corresponding to the aforementioned layers. For the most part, the trends seen in Figure 3 are reflected in Figure 4, with more warmer temperatures in the model at the 0–150 m layer, except for a string of cool temperatures roughly around the Northern Balearic front, just south of the typical convection patch in the wintertime. Below 150 m, there is less consistent structure in the bias, and the overall bias is reduced. The salinity bias is generally less spatially structured on all three levels, however the model is generally fresher than observations, with some saline pockets around the North Balearic front.

The area bounded by the black box in Figure 4 covers the area from 42 to 42.5°N and 4.25 to 5°E. Later, in Section 4, the vertical column of water within this bounding box is spatially averaged to study the temporal trends of the model, and is therefore a relevant area to investigate for major biases. Within this box, the bias follows the trends found in Figure 3: fresher and warmer water at the surface and fairly accurate at the lower layers. As water mass characteristics are typically separated on the order of 0.1 kg/m³ or more (Quentin-Boris et al., 2020), an order of magnitude larger than our density biases, and as we look at the whole vertical column for our study, we therefore believe the model results to be representative enough for our purpose of studying deep convection over multiple years.

## 3.2. Mixed Layer Depth

When comparing the model to the Argo calculated MLDs (same procedure as for the temperature and salinity profiles), the model overestimates the mixed layer depth. Over the domain in Figure 5, the average positive bias is 25 and 44 m (therefore meaning overestimation of the MLD by the model) when the model MLD is compared to the observed, estimated MLD. Around the deep convection area, this bias can much larger (over 500 m), which is due to the timing of deep convection, as the MLD grows rapidly during the phenomena and can be out of phase when compared to the observations. This is confirmed by subplot (c) of Figure 5, which shows almost all of the monthly variability in the bias occurring during the months of January, February, and March, when large vertical mixing can occur. This suggests that the model is fairly representative of the MLD but naturally struggles with the timing of deep convection, as it is a chaotic phenomenon (Waldman et al., 2018).

## 3.3. Wind

Figure 6 shows the bias of the WRF-ORCHIDEE atmospheric forcing data set used to force our ocean simulations and the QuikSCAT observed surface wind speeds. The WRF-ORCHIDEE data was daily averaged to be compared one-to-one with the QuikSCAT data. In general, in the northwest Mediterranean, the atmospheric forcing data underestimates the wind speeds by an average of -0.37 m/s, which is visible in Figure 6a, with an average RSME of 1.48 m/s. The same is seen in the differences between the 95th percentile of the forcing data and QuikSCAT, with an average difference of -1.13 m/s. Overall, the forcing data underestimates the wind at every percentile (shown in Figure S4 in Supporting Information S1), but not to an incredible extent. For the domain shown, the average correlation coefficient is 0.85, showing a high temporal correlation between the forcing data and QuikSCAT. As will be seen later in Section 4 and in Figure 11, the forcing data is reasonably representative as differences over 3–4 m/s in wind speed are needed to have a significant impact on the affected surface fluxes.

# 4. Results and Discussion

## 4.1. Seasonal Contribution

The center of the minimum point of stratification in the Gulf of Lion varies from year to year, which can be seen in Figure 5b. To compare the destratification from different years to each other, an area roughly encompassing the minima while maintaining decent coverage of the deepening MLD (again see Figure 5b) over the 20 years was averaged: a box with the limits of 42–42.5°N and 4.25–5°E (this box was first mentioned in Section 3.1). A video provided in the Supporting Information (Movie S1) shows the majority of the destratification in the winter time occurring in the zone. Seven years featured deep convection events in the model results of the control set: 1999, 2000, 2005, 2009, 2011, 2012, and 2013, shown by the significantly deep mixed layer depths (MLD) (in Figure 7; the years are highlighted with green text). This is in agreement with Somot et al. (2016) but in disagreement with observations shown in Bosse et al. (2021) and Houpert et al. (2016). Observations showed deep convection also occurred in 2010, but, as seen in Figures 8 and 10, our results show similar levels of stratification for 2010 as the
<figure><img src="image_8.png" /><figcaption>Figure 7.</figcaption></figure>
adjacent years, therefore capturing some of the behavior despite deep convection not occurring in the model. Our years of deep convection had the lowest stratification levels during convection, according to the stratification index (Figure 8), as expected (for the rest of the article, deep-convection years will refer to the deep-convection years found in the model results). Ideally, the stratification would be zero to denote a deep convection event. However, due to the area-averaging, some still stratified columns or laterally advected stratified waters are captured, resulting in some remaining stratification at the SI minimum for years with deep convection (hence why we tried to keep the bounding box smaller). This is particularly apparent for the year of 2009, a deep-convection year, that has some remaining stratification larger than the following years, due to 2009 having a deep-convection zone with a relatively small horizontal extent (see Movie S1).

The lack of deep convection in the seasonal set of simulations is immediately noticeable; the MLD for the seasonal runs never exceeded 200 m (Figure 7), except for 2005 where it reached a maximum of intermediate mixing at 586 m (spatially averaged; it reaches below 1,000 m for the 95th percentile of the region. This is shown in Figure S5 in Supporting Information S1). This confirms that the Mistral component is necessary for deep convection, as found for the winter of 2013 in Keller et al. (2022). However, there is a large variability of the seasonal stratification, $SI_S$. For example, for the winter of 2000 (referring to the winter spanning 1999–2000), the seasonal stratification closely follows the total stratification ($SI_S + \delta SI$), whereas the next winter, the winter of 2001, the seasonal stratification diverges quite strongly in February 2001 and remains diverged until June 2001 (Figure 8). To compare the variability between the different years, the seasonal and Mistral contributions, $\Omega SI_S$ and $\Omega \delta SI$, respectively, are determined according to Figure 9. The contributions are determined at the time where the total stratification reaches a minimum, $t_{SI_{min}}$, as this is where deep convection occurs in the years that feature a deep-convection event. This allows us to separate the contribution to destratification of each timescale:

$$ \Omega SI_S = SI_{S,max} - SI_S(t = t_{SI_{min}}) \qquad (12) $$

$$ \Omega \delta SI = SI_S(t = t_{SI_{min}}) - SI_{min} \qquad (13) $$

The maximum $SI_S$, noted as $SI_{S, max}$, is used as the reference point for the maximum stratification, at $t_{SI_{S,max}}$, as the seasonal stratification maximum is the overall stratification that both the Mistral and seasonal atmospheric change must overcome to cause deep convection. Consequently, the time $t_{SI_{S,max}}$ is taken to be the time the preconditioning phase begins, and the time $t_{SI_{min}}$ where it ends.

The varying levels of contributions and maximum levels of seasonal stratification are displayed in Figure 10 (with values given in Table S10 in Supporting Information S1). We can see in Figure 10a that the years with deep convection have $SI_{S, max}$ levels that are below average for the 20 year period ($1.69 m^2 s^{-2}$ for deep-convection (DC) years vs. $1.81 m^2 s^{-2}$ for all years, as seen in Table 3, first row, first column; DC years are denoted by the hatching in Figure 10). If we look at the separated contributions in subplots (b) and (c) of the same figure, the years with deep convection typically feature higher than average levels of destratification coming from the seasonal contribution, $\Omega SI_S$, with most of the destratification in 2012 coming from the seasonal timescale. For the Mistral timescale contributions, $\Omega \delta SI$, (see Figure 10d) years with deep convection also saw above average levels, except for the year of 2012. A key note of interest is the average levels of contribution from the two timescales. On average for all years, the seasonal timescale provides more of the required destratification necessary for deep convection than the Mistral timescale: 46% versus 28% (Table 3 second row, first column). This agrees with the results of Keller et al. (2022). If we compare their relative contribution to destratification in terms of their sum ($\Omega SI_S + \Omega \delta SI$), shown in Figure 10e and provided in the third row of Table 3, $\Omega SI_S$ provides more than half of the destratification (average of 61% for all years). When separated into DC and non-deep-convection (NDC) years, the averages remain above half, at 68% and 58%, respectively. This again confirms the importance of the seasonal destratification component, regardless of the occurrence of deep convection.

**Figure 7.** Averaged mixed layer depth of the area marked by the black box in Figure 4 for the 20 years, calculated by taking the depth to which the diffusivity is less than $5 \times 10^{-4} m^2/s$. The averaged mixed layer is fairly representative of the area, as seen in Figure S5 in Supporting Information S1, which shows the 95th percentile and maximum MLD for this area as well. In Figure S5 in Supporting Information S1, the mean, 95th percentile, and max show similar values, especially during deep convection years. The red circle labels the first date at which the MLD is deeper than 250 m and the red triangle marks the first main maximum depth for the deep-convection years. Mistral events are shown with the colored green shading.
<figure><img src="image_9.png" /><figcaption>**Figure 8.** The stratification index averaged over the area marked by the black box in Figure 4 for the 20 years, with the control run, *SI*<sub>S</sub> + *δSI*, in black and the seasonal run, *SI*<sub>S</sub>, in blue. The difference between the control and seasonal stratification index, *δSI*, is shown with a dashed red line with a separate scale. Mistral events are shown with the colored green shading.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Figure 9.** The stratification index for the winter of 1994 to demonstrate how the contributions from the different timescales are calculated.</figcaption></figure>

Analyzing the quantities in Table 3 a little further, we see that the average seasonal contribution, $\Omega SI_s$, is greater for DC years alone than all of the years combined or NDC years. This is true for the raw contributions (DC vs. overall and NDC years: 1.05 vs. 0.83 and 0.71 $m^2s^{-2}$) and relative contributions (62% vs. 46% and 37% relative to $SI_s$, max and 68% vs. 61% and 57% relative to $\Omega SI_s + \Omega \delta SI$; all quantities found in Table 3 first column). Compounding this with the lower than average $SI_s$, max for DC years as mentioned above, means that DC years have more destratification from the seasonal timescale and less stratification overall to overcome for deep convection to occur.

Conversely, the distinction between DC, NDC, and the overall average for the Mistral contributions, $\Omega \delta SI$, is less clear, neither in relative nor absolute terms. The only distinction is found in the relative contribution (with regard to $\Omega SI_s + \Omega \delta SI$) of the high frequency component of the Mistral, which is less

during DC than NDC years (32% vs. 43%; Figure 10e and the bottom row of Table 3). The Mistral contribution is also more variable during NDC years than DC years in relative and absolute terms (Table 3, second column). However, for the raw and normalized contributions relative to $SI_s$, max, there is little difference between DC and NDC years (0.48 vs. 0.50 $m^2s^{-2}$ – 2% and 29% and 28%, respectively, found in Table 3, first and second row, first column). Despite this variability, the high frequency component of the Mistral doesn't vary much between the various groupings of years (DC, NDC, and all years).

### 4.1.1. Components of the Seasonal Contribution

As the variability of the seasonal contribution, $\Omega SI_s$, plays a key role in the occurrence of deep convection as seen above, it was separated into the different surface heat flux components, as described by Equation 9, with $t_0 = t_{SI_s, max}$ and $t_1 = t_{SI_s, min}$. The distributions of the different flux components over the years are shown in Figure 11a, with DC years colored in blue and NDC years colored in red. What Figure 11a conveys, is that the years with increased latent, $Q_E$, and sensible, $Q_H$, heat fluxes during the preconditioning period are the years with deep convection. This is seen by the differences in the mean values for each subgroup ($\overline{DC} - \overline{NDC}$): 0.04, -0.11, -0.17, and -0.43 $m^2s^{-2}$ for $Q_{SW}$, $Q_{LW}$, $Q_H$, and $Q_E$, respectively. Initially, it appears that the longwave upward radiation also acts as an indicator for years with deep convection. However, if we normalize these differences by the average value of all of the years for each component ($\overline{(\overline{DC} - \overline{NDC})}/\overline{All}$), then we can easily identify the sensible and latent heat fluxes as the main indicators: 0.03, 0.08, **0.43**, and **0.33** (in the same order as the previous list). The resulting p-values of 0.85, 0.40, 0.00, 0.01 between the DC and NDC years support this conclusion as well, with only the sensible and latent heat flux rejecting the null hypothesis of both DC and NDC years coming from the same distribution.

To determine which atmospheric component drove the differences evident in the latent and sensible heat fluxes, $\Delta q$, $\Delta \theta$, and the wind speed ($|\Delta \vec{u}| \approx |\vec{u}_z|$, as the sea surface current is typically small relative to the wind speed; typically $\mathcal{O}(0.1-1\text{ m/s})$ versus. $\mathcal{O}(10\text{ m/s})$, respectively) were ensemble averaged for DC and NDC years (blue and red, respectively). These values were selected as they are the atmospheric components found in Equation 4 used to calculate the latent and sensible heat fluxes. The ensemble averaging is shown in Figure 11 subplots (b), (c), and (d), respectively. While there are differences in both $\Delta q$ and $\Delta \theta$ between DC and NDC years, the wind speed, $|\vec{u}_z|$ is the main differentiator between the two groups of years. This is quantitatively seen with daily t tests, comparing the above variables between DC and NDC years. The hatched areas in Figures 11a–11c denote where the daily p-values are below 5%, indicating where the null hypothesis, assuming the DC and NDC variables are from the same population, fails; these areas are larger and more frequent in the wind speed than in the other variables.

To better demonstrate $|\vec{u}_z|$ as the main differentiator, a sensitivity analysis was performed by estimating the change in destratification due to the latent and sensible heat fluxes, with either the DC or NDC ensemble averaged values for $\Delta q$, $\Delta \theta$, and $|\vec{u}_z|$. Using Equation 9, the change in the estimated destratification, denoted by subscript “Est,” due to changes in $Q_E$ and $Q_H$, can be calculated as:
<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

**Figure 10.** The seasonal maximum stratification and minimum control stratification are shown in subplot (a). The seasonal and Mistral contributions are shown in (b) and (c) (normalized in (c)). (d) shows the normalized Mistral contribution to destratification. (e) shows the contributions normalized to their sum. The solid lines are the means of the indicated variables. The dashed lines are the trend lines (none break the null hypothesis of a zero slope: *p*-value *≈* 0.05). The overlaid hatching signifies years that experienced deep convection.
<figure><img src="image_16.png" /><figcaption>**Figure 11.** Distributions of the different flux components making the stratification change in SI<sub>s</sub>, determined using Equation 9 per component, i, are shown in subplot (a). DC and NDC stand for deep-convection and non-deep-convection, respectively. Subplots (b)–(d) show the ensemble averaged (discarding February 29th from leap years) driving components of the flux bulk formulae in Equation 4, Δq, Δθ, and |u̇<sub>z</sub>|. In subplot (a), the solid lines are the mean of the indicated group; in subplots (b)–(d), the solid lines are the ensemble mean and the shading is plus/minus one standard deviation. The hatching denotes the days that break the null hypothesis that the DC and NDC variables are from the same population, with p-values under 5%.</figcaption></figure>

$$ \Delta SI_{Est,i,j} = 10^{-9} \times \int_{t_0}^{t_1} \Delta Q_{i,j} dt \quad (14) $$

Where ΔQ<sub>i,j</sub> = Q<sub>i,j</sub> - Q<sub>i,Ref</sub>, i is either E or H for the latent and sensible heat flux, respectively, and j is either Δq, Δθ, or |u̇<sub>z</sub>|. Here, j stands for the variable changed to the DC ensemble averaged value (denoted by the subscript DC), setting the remaining variables to the NDC ensemble averaged values (denoted by the subscript NDC). Q<sub>i,Ref</sub> has all variables set to the NDC ensemble averaged values. For example, ΔQ<sub>E,Δq</sub> would be:
<table><caption>**Table 3**

Statistics for the Non-Normalized and Normalized $\Omega SI_s$ and $\Omega \delta SI$
Contributions From Figures 10b–10e</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th colspan="3">Mean</th><th colspan="3">Std</th><th colspan="3">Max</th></tr><tr><th>ΩSI<sub>s</sub></th><th>ΩδSI</th><th>SI<sub>s, max</sub></th><th>ΩSI<sub>s</sub></th><th>ΩδSI</th><th>SI<sub>s, max</sub></th><th>ΩSI<sub>s</sub></th><th>ΩδSI</th><th>SI<sub>s, max</sub></th></tr></thead><tbody><tr><td colspan="10">Non-normalized m<sup>2</sup>s<sup>-2</sup></td></tr><tr><td>DC</td><td>1.05</td><td>0.48</td><td>1.69</td><td>0.29</td><td>0.09</td><td>0.19</td><td>1.69</td><td>0.59</td><td>2.10</td></tr><tr><td>NDC</td><td>0.71</td><td>0.50</td><td>1.88</td><td>0.34</td><td>0.21</td><td>0.41</td><td>1.55</td><td>0.97</td><td>3.00</td></tr><tr><td>All</td><td>0.83</td><td>0.49</td><td>1.81</td><td>0.36</td><td>0.18</td><td>0.36</td><td>1.69</td><td>0.97</td><td>3.00</td></tr><tr><td colspan="10">Normalized · (/SI<sub>s, max</sub>)</td></tr><tr><td>DC</td><td>0.62</td><td>0.29</td><td>—</td><td>0.12</td><td>0.06</td><td>—</td><td>0.81</td><td>0.35</td><td>—</td></tr><tr><td>NDC</td><td>0.37</td><td>0.28</td><td>—</td><td>0.12</td><td>0.11</td><td>—</td><td>0.52</td><td>0.52</td><td>—</td></tr><tr><td>All</td><td>0.46</td><td>0.28</td><td>—</td><td>0.17</td><td>0.10</td><td>—</td><td>0.81</td><td>0.52</td><td>—</td></tr><tr><td colspan="10">Normalized ( · / (ΩSI<sub>s</sub> + ΩδSI))</td></tr><tr><td>DC</td><td>0.68</td><td>0.32</td><td>—</td><td>0.08</td><td>0.08</td><td>—</td><td>0.59</td><td>0.15</td><td>—</td></tr><tr><td>NDC</td><td>0.57</td><td>0.43</td><td>—</td><td>0.17</td><td>0.17</td><td>—</td><td>0.18</td><td>0.22</td><td>—</td></tr><tr><td>All</td><td>0.61</td><td>0.39</td><td>—</td><td>0.16</td><td>0.16</td><td>—</td><td>0.18</td><td>0.15</td><td>—</td></tr></tbody></table>

Note. DC stands for deep convection and NDC stands for non-deep convection.

$$ \begin{aligned} \Delta Q_{E, \Delta q} &= Q_{E, \Delta q} - Q_{E, Ref} \\ &= \rho_{a,0} \Lambda C_E (\Delta q_{DC}) |\vec{u}_z|_{NDC} - \rho_{a,0} \Lambda C_E (\Delta q_{NDC}) |\vec{u}_z|_{NDC} \end{aligned} \quad (15) $$

We can then determine the direct influence that the DC ensemble averaged values for $\Delta q$, $\Delta \theta$, and $|\vec{u}_z|$ have on destratification. The results of this analysis are found in Table 4. As $|\vec{u}_z|$ influences both $Q_E$ and $Q_H$, it easily makes a larger difference in terms of destratification than either $\Delta q$ or $\Delta \theta$: $-0.42 \text{ m}^2\text{s}^{-2}$ versus $-0.066 \text{ m}^2\text{s}^{-2}$ and $-0.042 \text{ m}^2\text{s}^{-2}$, respectively.

The source of this difference in wind speed between DC and NDC years obfuscates the distinction between seasonal and Mistral contributions, however. The filtering, as discussed in Keller et al. (2022) and in Section 2.2, primarily removes the high frequency component of the Mistral (resulting filtered heat fluxes are shown in Figure S2 in Supporting Information S1). However, as also pointed out, the Mistral has a low frequency seasonal component as well, with more frequent and stronger Mistrals occurring in winter versus summer (see Givon et al. (2021) for a more complete analysis). With the moving average window, this low frequency component is partially filtered out, removing some of the Mistral's low frequency component (when viewed in the spectral domain), however part of it still remains. This remaining low frequency component of the Mistral is the overall increase in the mean wind speed during the winter months due to more frequent Mistral

events and imprints itself in the seasonal forcing. The percentage of the preconditioning days ($t_{SI_{s,max}}$ to $t_{SI_{min}}$) that feature a Mistral event supports this observation, with DC years at 34% and NDC years at 29%.

## 4.2. Discussion

As mentioned in the introduction, in Keller et al. (2022) a few questions were posed that couldn't be answered by the scope of that study. We will readdress them here.

### 4.2.1. Does the Mistral Trigger Deep Convection, or Does the Seasonal Change Trigger It?

To determine if the Mistral or seasonal change triggered deep convection in our study, or which triggers the intense vertical mixing, we first located the main growth phase of the MLD during deep convection. This was done to focus the attention specifically on the time period where vertical mixing occurs, rather than over the entire preconditioning period (from $SI_{s, max}$ to $SI_{min}$). The main growth phase was chosen to be the first point in time at which the MLD became deeper than 250 m (labeled by a red circle in Figure 7 for DC years) to the point at which the MLD reaches its first maximum (first if two major peaks were present, such as for the years of 2009, 2012, and 2013, otherwise the overall maximum was used; labeled by a red triangle in Figure 7 for the same years). Then the ratio of the averaged gradient, with respect to time, of $\delta SI$ and $SI_S$ ($\partial_t \delta SI / \partial_t SI_S$) was

computed for this growth phase for each DC year (visualized in Figure S6 for 1999 in Supporting Information S1). By comparing these two gradients of stratification, we can directly determine which surface forcing, the Mistral or seasonal, contributes more to the destabilization of the vertical column during the intense vertical mixing, as the temporal gradient is proportional to the surface net heat flux (effectively generalizing Equation 8). This thereby determines which one is driving this vertical mixing. The results are shown in Figure 12. The years of 2000, 2009, and 2013 saw a larger destratifying contribution from the high frequency part of the Mistral than the seasonal forcing, with ratios greater than unity: 1.45, 4.71, and 2.15, respectively. This demonstrates that the Mistral was the main triggering component for these years. However, for the years of 1999, 2005, 2011, and 2012, the seasonal component was the main triggering agent, with ratios less than unity: 0.43, 0.40, 0.18, and 0.05, respectively. This means

<table><caption>**Table 4**

Estimated Changes in Destratification Due To Changing One Variable at a Time (Between $\Delta q$, $\Delta \theta$, and $|\vec{u}_z|$) to DC Versus NDC Ensemble Averaged Values, Utilizing Equation 14</caption><thead><tr><th>j</th><th>ΔSI<sub>Est,E,j</sub> m<sup>2</sup>/s<sup>2</sup></th><th>ΔSI<sub>Est,H,j</sub> m<sup>2</sup>/s<sup>2</sup></th><th>Total m<sup>2</sup>/s<sup>2</sup></th></tr></thead><tbody><tr><td>Δq</td><td>-0.066</td><td>-0</td><td>-0.066</td></tr><tr><td>Δθ</td><td>-0</td><td>-0.042</td><td>-0.042</td></tr><tr><td>|u̇<sub>z</sub>|</td><td>-0.32</td><td>-0.10</td><td>-0.42</td></tr></tbody></table>

Note. Note, the saturation humidity is based on sea surface temperature, which means keeping the temperature at NDC ensemble averaged values is technically non-physical, as the saturation humidity would change with a different air temperature.
<figure><img src="image_17.png" /><figcaption>**Figure 12.** The averaged gradients, with respect to time, of $\delta SI$ and $SI_s$ during the main growth phase of vertical mixing for each DC year: $\partial_t \delta SI$ and $\partial_t SI_s$ as the green and blue bars, respectively (the overbars have been omitted in the figure legend for readability). The ratio of the former over the latter is shown atop the graph. The bars' values are also displayed at the bottom of the corresponding bars.</figcaption></figure>

both the Mistral and seasonal component trigger deep convection in roughly equal amounts of our studied DC years, and neither dominate.

### 4.2.2. Does the Maximum Seasonal Stratification Play a Role in Deep Convection?

According to our results, the maximum $SI_s$ does play a role. As previously pointed out, DC years are almost entirely years with a lower than average $SI_{s,max}$ (except for 2012). Which is intuitive, as a larger maximum of $SI_s$ means that both the seasonal component and Mistral component must overcome a larger amount of stratification to form deep convection. However, more importantly, years with above average $\Omega SI_s$ are more often than not, DC years. We saw the origin of this difference in Figure 11 and Table 3, in the difference of wind speed. This means that the seasonal contribution to destratification, through the wind speed, has a particularly important role in the overall destratification of the Gulf of Lion, as well as the $SI_{s,max}$ it must overcome.

A similar question could also be asked: does the previous year's level of stratification affect the proceeding year? Five of the seven DC years occurred adjacently in our simulations: the years of 1999 and 2000 occurred together and the years 2011–2013 occurred together as well. Otherwise, the two remaining years were in between two NDC years. This seems to suggest DC years occur consecutively, which intuitively makes sense, as the water column following a deep convection event will have had a significant amount of heat removed from it (resulting in buoyancy loss, driving destratification). This heat must be re-injected into the water column to restratify it, whereas years with persisting stratification don't need this initial addition of heat. The converse is seen in Margirier et al. (2020), where deep convection is inhibited after 2013 and leads to a warming and more saline intermediate water in the gulf, preventing deep convection in these years.

Returning to our results, however, in terms of stratification (through the stratification index and contributions derived from the $SI$), there doesn't appear to be any discernible pattern or trend for the 20 year period we studied (see Figure 10), as the trend in Margirier et al. (2020) occurs after 2013 (Figure 5c of their work). A larger scoped study that investigates additional features, such as the composition of the formed dense water masses (e.g., the saltier dense water formed during the 2005 deep convection event (Herrmann et al., 2010)), the long term trends of said composition (Houpert et al., 2016), or the changes in the Mediterranean Sea circulation (Amitai et al., 2021), may be able to provide more answers. For example, the study of Parras-Berrocal et al. (2022) found that increasingly saline Levantine Intermediate Water and freshening Inflow Atlantic Water at the Strait of Gibraltar lead to increasing stratification in the Gulf of Lion for climatic scenario runs til the year of 2100.
<figure><img src="image_18.png" /><figcaption>**Figure 13.** Distributions of the destratification incurred by Mistral events during DC years. $\Delta\delta SI$ is calculated using Equation 10 per event k. The solid vertical lines are the means of the indicated groups.</figcaption></figure>

### 4.2.3. Do Any of the Deep Convection Events Occur After the Minimum in the Seasonal Stratification?

The third question, broken down into a few separate yet related questions, poses: does the timing of the $SI_S$ minimum matter? Can the Mistral, $\delta SI$, overcome the restratifying $SI_S$? Or, in other words, do any of the deep convection events occur after the $SI_S$ minimum?

For our results, three of the seven DC years (2009, 2012, and 2013) experienced a control $SI$ minimum that occurred after the $SI_S$ minimum (vertical dashed lines in Figure 7). In each of these 3 years, according to the MLD (Figure 7), deep convection ceased temporarily between the control and seasonal stratification minimum. Then deep convection resumed with an additional peak in the MLD before the control $SI$ reached it's minimum. This means that the Mistral can overcome a restratifying $SI_S$ to continue an already triggered deep convection. This was observed in 2010 by Houpert et al. (2016) at the end of February and beginning of March, when the ocean column briefly restratified with warming surface fluxes before cooling fluxes returned and resumed vertical mixing. However, it is unclear whether or not the Mistral can trigger deep convection after the seasonal minimum, as our model results don't feature such an example.

While a larger data set of deep convection events will be required to more definitively answer this question, we can infer that the case of triggering deep convection after or continuing it beyond the $SI_S$ minimum will be rarer than the case of the control $SI$ minimum occurring before the $SI_S$ minimum. This is due to a weakening contribution from Mistral events as the preconditioning period occurs. Equation 19 of Keller et al. (2022) shows that succeeding Mistral events need to be stronger than the current level of destratification to cause more destratification. This is compounded by the fact that the average wind speed for Mistral days drop after February (not shown), additionally reducing its strength. When the year transitions out of the preconditioning period into the summer (essentially after the $SI_S$ minimum), Mistral events destratify less because the water column has already incurred a significant amount of destratification. We can see this change by looking at the destratification caused by individual Mistral events depending on their timing with Equation 10. The results are shown in Figure 13. Events before the $SI_S$ minimum exhibit a wider spread in terms of destratification, but also have a mean destratification ($-0.021 \text{ m}^2\text{s}^{-2}$) that is less than the events that occur after the minimum ($-0.001 \text{ m}^2\text{s}^{-2}$). This limits the likelihood that Mistral events can overcome a restratifying $SI_S$.

## 5. Conclusions

Our study investigated deep convection in the Gulf of Lion over a 20 year period, using the NEMO ocean model, in the NEMOMED12 configuration, forced by filtered and unfiltered RegIPSL WRF/ORCHIDEE atmospheric data. By looking at the difference between the two sets of ocean simulation results forced by the two different forcings, we could extract the effect the Mistral and seasonal atmospheric change had on the annual stratification cycle of the Gulf of Lion. The control model results represented reality fairly well with respect to Argo and CTD
profiling, with the same being concluded for the comparison between the control model MLDs and the Argo MLDs. While deep convection occurs in only seven of the 20 years in the model results, whereas it occurs in eight of the 20 years in observations (Bosse et al., 2021; Houpert et al., 2016), we were able to extract information regarding the impact of the seasonal atmospheric change on destratification.

We found the seasonal contribution of the atmospheric forcing to be the main driver in terms of destratification during the preconditioning period, with it being larger during DC years. In addition to this, DC years featured lower levels of maximum stratification to overcome for deep convection to occur, requiring less forcing needed to destratify the water column. The high frequency Mistral component, however, did not vary much between DC and NDC years but was determined to still be necessary for deep convection to occur.

When breaking down what causes destratification in the seasonal contribution, we found the latent and sensible heat fluxes to be most important components, shifting more negative during DC years. It was then found that the differences in the latent and sensible heat fluxes between DC and NDC years were caused by increased wind speeds during DC years. These increased wind speeds themselves were caused by the low frequency aspect of Mistral winds, which is imprinted in the seasonal forcing, with more events occurring during the winters with deep convection.

When addressing the questions asked by Keller et al. (2022), we found that the Mistral and seasonal atmospheric change roughly trigger deep convection an equal number of times. It was also determined that the maximum seasonal stratification is an important quantity as it is the amount of stratification the seasonal and Mistral contributions must overcome to cause deep convection. Additionally, the Mistral contribution can overcome a restratifying seasonal stratification to extend deep convection, however it is unlikely it can trigger deep convection after the seasonal stratification minimum.

Our study shows the importance of the seasonal atmospheric change and its drivers on the deep convection cycle of the Gulf of Lion. Studying the effects of climate change on the seasonal atmospheric forcing and vertical composition of the Gulf of Lion waters will be necessary to understand the evolution of this phenomena, highlighting the sensitivity of the Mediterranean to climate change (Adloff et al., 2015; Herrmann et al., 2008, 2014; Parras-Berrocal et al., 2022; Somot et al., 2006).

## Appendix A: Determining $\delta F_k$

The relationship between $\delta SI$ and individual Mistral events, $k$, can be described as the following (from Keller et al. (2022) Equation A35):

$$ \delta SI_k(t) = \begin{cases} \left[ \delta SI_{k-1}(t_k) + \frac{D^2}{2} \frac{\delta F_k}{\alpha_d} (1 - e^{\alpha_d(t-t_k)}) \right] e^{-\alpha_d(t-t_k)} & [t_k, t_k + \Delta t_k) \\ \left[ \delta SI_{k-1}(t_k) + \frac{D^2}{2} \frac{\delta F_k}{\alpha_d} (1 - e^{\alpha_d \Delta t_k}) \right] e^{(\alpha_a - \alpha_d)\Delta t_k - \alpha_a(t-t_k)} & [t_k + \Delta t_k, t_k + \Delta \tau_k) \end{cases} \quad (A1) $$

With the restoring coefficients, $\alpha_a$ and $\alpha_d$, the duration, $\Delta t_k$, and the period, $\Delta \tau_k$, the strength of each Mistral event, $\delta F_k$, can be determined. If we take Equation A1 and note the value of $\delta SI_{k-1}$ to be the same as $\delta SI_k(t_k)$ at the beginning of an event, we can simplify the equation in the following steps:

$$ \delta SI_k(t_k) = \delta SI_{k-1} \quad (A2) $$

$$ \delta SI_k(t_k + \Delta t_k) = \delta SI_{k-1}e^{-\alpha_d \Delta t_k} + \frac{D^2}{2} \frac{\delta F_k}{\alpha_d}(e^{-\alpha_d \Delta t_k} - 1) \quad (A3) $$

And then solve for $\delta F_k$:

$$ \delta F_k = \frac{2(\delta SI_k(t_k + \Delta t_k) - \delta SI_k(t_k) e^{-\alpha_d \Delta t_k})\alpha_d}{(e^{-\alpha_d \Delta t_k} - 1)D^2} \quad (A4) $$
Tables S1 through S9 in Supporting Information S1 show the resulting values for each Mistral event over the 20 year period, calculated based on the δSI from our modeling results.

## Data Availability Statement

The NEMO simulation data used in this article is available at https://doi.org/10.57745/D7ACY9; last accessed 29 March 2024 (Keller, 2024). The Mistral event data is available at https://medcyclones.eu/data/; last accessed 23 August 2022 (Givon, 2021; Givon et al., 2021; Keller et al., 2022). The RegIPSL WRF/ORCHIDEE atmospheric forcing data is available from the authors of Guion et al. (2021b) upon request (Guion et al., 2021a). The Argo and CTD data is available through the Coriolis database, available at https://www.coriolis.eu.org/Data-Products/Data-selection; last accessed: 23 August 2022 (Argo and CTD vertical profiles, 2022). The Argo mixed layer depth data is available at https://mixedlayer.ucsd.edu; last accessed: 22 February 2023 (Holte et al., 2017a, 2017b). The QuikSCAT is presented at https://cersat.ifremer.fr/Data/Latest-products/Reprocessing-Daily-QuikSCAT-Surface-Wind-Fields; last accessed: 22 February 2023, however at the time of writing, availability is limited due the migration of its data (see https://podaac.jpl.nasa.gov/cloud-datasets/migration; last accessed 16 June 2023) (NASA, 2023). The RegIPSL model can be found at https://gitlab.in2p3.fr/ipsl/lmd/intro/regipsl/regipsl; last accessed: 26 August 2022 (InTro LMD, IPSL, 2022). The NEMO model can be found at https://www.nemo-ocean.eu/; last accessed: 16 August 2022 (Madec et al., 2022). The code used to perform the analysis and produce the plots are available at https://gitlab.com/dkllrjr/jgr_oceans_untanglingolumconv_20yrs_code; last accessed: 24 August 2022 (Keller, 2023).

## Acknowledgments

This work is a contribution to the HyMeX program (Hydrological cycle in the Mediterranean EXperiment) through INSU-MISTRALS support and the MED-CORDEX program (Coordinated Regional climate Downscaling EXperiment—Mediterranean region). It was also supported by a joint CNRS—Weizmann Institute of Science collaborative project.

## References

Adloff, F., Somot, S., Sevault, F., Jordà, G., Aznar, R., Déqué, M., et al. (2015). Mediterranean Sea response to climate change in an ensemble of twenty first century scenarios. *Climate Dynamics*, *45*(9–10), 2775–2802. https://doi.org/10.1007/s00382-015-2507-3

Amitai, Y., Ashkenazy, Y., & Gildor, H. (2021). The effect of the source of deep water in the Eastern Mediterranean on Western Mediterranean intermediate and deep water. *Frontiers in Marine Science*, *7*, 615975. https://doi.org/10.3389/fmars.2020.615975

Argo and CTD vertical profiles. (2022). Argo and CTD vertical profiles [Dataset]. Retrieved from https://www.coriolis.eu.org/Data-Products/Data-selection

Balmaseda, M. A., Trenberth, K. E., & Källén, E. (2013). Distinctive climate signals in reanalysis of global ocean heat content. *Geophysical Research Letters*, *40*(9), 1754–1759. https://doi.org/10.1002/grl.50382

Beuvier, J., Béranger, K., Brossier, C. L., Somot, S., Sevault, F., Drillet, Y., et al. (2012). Spreading of the western Mediterranean deep water after winter 2005: Time scales and deep cyclone transport. *Journal of Geophysical Research*, 117(C7). https://doi.org/10.1029/2011JC007679

Bosse, A., Testor, P., Damien, P., Estournel, C., Marsaleix, P., Mortier, L., et al. (2021). Wind-forced submesoscale symmetric instability around deep convection in the northwestern Mediterranean Sea. *Fluid*, *6*(3), 123. https://doi.org/10.3390/fluids6030123

Bosse, A., Testor, P., Houpert, L., Damien, P., Prieur, L., Hayes, D., et al. (2016). Scales and dynamics of submesoscale coherent vortices formed by deep convection in the northwestern Mediterranean Sea. *Journal of Geophysical Research: Oceans*, *121*(10), 7716–7742. https://doi.org/10.1002/2016jc012144

Coppola, L., Prieur, L., Taupier-Letage, I., Estournel, C., Testor, P., Lefevre, D., et al. (2017). Observation of oxygen ventilation into deep waters through targeted deployment of multiple ARGO-O2FLOATS in the north-western Mediterranean Sea in 2013. *Journal of Geophysical Research: Oceans*, *122*(8), 6325–6341. https://doi.org/10.1002/2016jc012594

Drobinski, P., Alonzo, B., Basdevant, C., Cocquerez, P., Doerenbecher, A., Fourrié, N., & Nuret, M. (2017). Lagrangian dynamics of the mistral during the HyMeX SOP2. *Journal of Geophysical Research: Atmospheres*, *122*(3), 1387–1402. https://doi.org/10.1002/2016jd025530

Drobinski, P., Ducrocq, V., Alpert, P., Anagnostou, E., Béranger, K., Borga, M., et al. (2014). HyMeX: A 10-year multidisciplinary program on the Mediterranean water cycle. *Bulletin of the American Meteorological Society*, *95*(7), 1063–1082. https://doi.org/10.1175/bams-d-12-00242.1

Estournel, C., Testor, P., Damien, P., D'Ortenzio, F., Marsaleix, P., Conan, P., et al. (2016). High resolution modeling of dense water formation in the north-western Mediterranean during winter 2012–2013: Processes and budget. *Journal of Geophysical Research: Oceans*, *121*(7), 5367–5392. https://doi.org/10.1002/2016jc011935

Firing, E., Filipe, Barna, A., & Abernathey, R. (2021). Teos-10/gsw-python:v3.4.1.post0. Zenodo. https://doi.org/10.5281/ZENODO.5214122

Flamant, C. (2003). Alpine lee cyclogenesis influence on air-sea heat exchanges and marine atmospheric boundary layer thermodynamics over the western Mediterranean during a tramontane/mistral event. *Journal of Geophysical Research*, *108*(C2). https://doi.org/10.1029/2001jc001040

Frajka-Williams, E., Rhines, P. B., & Eriksen, C. C. (2014). Horizontal stratification during deep convection in the Labrador Sea. *Journal of Physical Oceanography*, *44*(1), 220–228. https://doi.org/10.1175/jpo-d-13-069.1

Givon, Y. (2021). Mistral data [Dataset]. A list of dates when mistral wind prevailed in the Gulf of Lion for a climatological period. https://medcyclones.eu/data/

Givon, Y., Keller, D., Jr., Silverman, V., Pennel, R., Drobinski, P., & Raveh-Rubin, S. (2021). Large-scale drivers of the mistral wind: Link to Rossby wave life cycles and seasonal variability. *Weather and Climate Dynamics*, *2*(3), 609–630. https://doi.org/10.5194/wcd-2-609-2021

Guion, A., Turquety, S., Polcher, J., Pennel, R., Bastin, S., & Arsouze, T. (2021a). Data from “Droughts and heatwaves in the Western Mediterranean: Impact on vegetation and wildfires using the coupled WRF-ORCHIDEE regional model (RegIPSL)” [Dataset]. *Climate Dynamics*, *58*(9–10), 2881–2903. https://doi.org/10.1007/s00382-021-05938-y

Guion, A., Turquety, S., Polcher, J., Pennel, R., Bastin, S., & Arsouze, T. (2021b). Droughts and heatwaves in the western Mediterranean: Impact on vegetation and wildfires using the coupled WRF-ORCHIDEE regional model (RegIPSL). *Climate Dynamics*, *58*(9–10), 2881–2903. https://doi.org/10.1007/s00382-021-05938-y
Hamon, M., Beuvier, J., Somot, S., Lellouche, J.-M., Greiner, E., Jordà, G., et al. (2016). Design and validation of MEDRYS, a Mediterranean Sea reanalysis over the period 1992–2013. *Ocean Science*, 12(2), 577–599. https://doi.org/10.5194/os-12-577-2016

Herrmann, M., Auger, P.-A., Ulses, C., & Estournel, C. (2017). Long-term monitoring of ocean deep convection using multisensors altimetry and ocean color satellite data. *Journal of Geophysical Research: Oceans*, 122(2), 1457–1475. https://doi.org/10.1002/2016jc011833

Herrmann, M., Estournel, C., Adloff, F., & Diaz, F. (2014). Impact of climate change on the northwestern Mediterranean Sea pelagic planktonic ecosystem and associated carbon cycle. *Journal of Geophysical Research: Oceans*, 119(9), 5815–5836. https://doi.org/10.1002/2014jc010016

Herrmann, M., Estournel, C., Déqué, M., Marsaleix, P., Sevault, F., & Somot, S. (2008). Dense water formation in the Gulf of Lions shelf: Impact of atmospheric interannual variability and climate change. *Continental Shelf Research*, 28(15), 2092–2112. https://doi.org/10.1016/j.csr.2008.03.003

Herrmann, M., Sevault, F., Beuvier, J., & Somot, S. (2010). What induced the exceptional 2005 convection event in the northwestern Mediterranean basin? Answers from a modeling study. *Journal of Geophysical Research*, 115(C12). https://doi.org/10.1029/2010jc006162

Herrmann, M., & Somot, S. (2008). Relevance of ERA40 dynamical downscaling for modeling deep convection in the Mediterranean Sea. *Geophysical Research Letters*, 35(4). https://doi.org/10.1029/2007g1032442

Holte, J., & Talley, L. (2009). A new algorithm for finding mixed layer depths with applications to Argo data and subAntarctic mode water formation. *Journal of Atmospheric and Oceanic Technology*, 26(9), 1920–1939. https://doi.org/10.1175/2009jtecho543.1

Holte, J., Talley, L. D., Gilson, J., & Roemmich, D. (2017a). An Argo mixed layer climatology and database. *Geophysical Research Letters*, 44(11), 5618–5626. https://doi.org/10.1002/2017gl073426

Holte, J., Talley, L. D., Gilson, J., & Roemmich, D. (2017b). An Argo mixed layer climatology and database [Dataset]. *Geophysical Research Letters*, 44(11), 5618–5626. https://doi.org/10.1002/2017gl073426

Houpert, L., de Madron, X. D., Testor, P., Bosse, A., D'Ortenzio, F., Bouin, M. N., et al. (2016). Observations of open-ocean deep convection in the northwestern Mediterranean Sea: Seasonal and interannual variability of mixing and deep water masses for the 2007–2013 period. *Journal of Geophysical Research: Oceans*, 121(11), 8139–8171. https://doi.org/10.1002/2016jc011857

Hurlburt, H., Brassington, G., Drillet, Y., Kamachi, M., Benkiran, M., Bourdallé-Badie, R., et al. (2009). High-resolution global and basin-scale ocean analyses and forecasts. *Oceanography*, 22(3), 110–127. https://doi.org/10.5670/oceanog.2009.70

InTro LMD, IPSL. (2022). The regional Earth system model of IPSL (RegIPSL) [Software]. Retrieved from https://gitlab.in2p3.fr/ipsl/lmd/intro/regipsl/regipsl

Keller, D. (2024). Data for untangling the mistral and seasonal atmospheric forcing driving deep convection in the Gulf of Lion: 1993–2013 [Dataset]. *Recherche Data Gouv*. https://doi.org/10.57745/D7ACY9

Keller, D., Jr. (2023). Repository of code used for this article [Software]. Retrieved from https://gitlab.com/dkllrjr/jgr_oceans_untanglingolumconv_20yrs_code

Keller, D., Jr., Givon, Y., Pennel, R., Raveh-Rubin, S., & Drobinski, P. (2022). Untangling the mistral and seasonal atmospheric forcing driving deep convection in the Gulf of Lion: 2012–2013. *Ocean Science*, 18(2), 483–510. https://doi.org/10.5194/os-18-483-2022

Krinner, G., Viovy, N., de Noblet-Ducoudré, N., Ogée, J., Polcher, J., Friedlingstein, P., et al. (2005). A dynamic global vegetation model for studies of the coupled atmosphere-biosphere system. *Global Biogeochemical Cycles*, 19(1). https://doi.org/10.1029/2003gb002199

Large, W. G., & Yeager, S. G. (2004). Diurnal to decadal global forcing for ocean and sea-ice models: The data sets and flux climatologies (Tech. Rep.). https://doi.org/10.5065/D6KK98Q6

Large, W. G., & Yeager, S. G. (2008). The global climatology of an interannually varying air–sea flux data set. *Climate Dynamics*, 33(2–3), 341–364. https://doi.org/10.1007/s00382-008-0441-3

Lebeaupin-Brossier, C., Béranger, K., Deltel, C., & Drobinski, P. (2011). The mediterranean response to different space–time resolution atmospheric forcings using perpetual mode sensitivity simulations. *Ocean Modelling*, 36(1–2), 1–25. https://doi.org/10.1016/j.ocemod.2010.10.008

Lebeaupin-Brossier, C., Béranger, K., & Drobinski, P. (2012). Sensitivity of the northwestern Mediterranean Sea coastal and thermohaline circulations simulated by the 1/12°-resolution ocean model NEMO-MED12 to the spatial and temporal resolution of atmospheric forcing. *Ocean Modelling*, 43–44, 94–107. https://doi.org/10.1016/j.ocemod.2011.12.007

Lebeaupin-Brossier, C., & Drobinski, P. (2009). Numerical high-resolution air-sea coupling over the Gulf of Lions during two tramontane/mistral events. *Journal of Geophysical Research*, 114(D10). https://doi.org/10.1029/2008jd011601

Ludwig, W., Dumont, E., Meybeck, M., & Heussner, S. (2009). River discharges of water and nutrients to the Mediterranean and black sea: Major drivers for ecosystem changes during past and future decades? *Progress in Oceanography*, 80(3–4), 199–217. https://doi.org/10.1016/j.pocean.2009.02.001

Madec, G. (2008). NEMO ocean engine, Note du Pôle de modélisation (Tech. Rep.). Institut Pierre-Simon Laplace.

Madec, G., Bourdallé-Badie, R., Chanut, J., Clementi, E., Coward, A., Ethé, C., et al. (2022). Nemo ocean engine [Software]. https://doi.org/10.5281/ZENODO.1464816

Madec, G., Chartier, M., & Crépon, M. (1991). The effect of thermohaline forcing variability on deep water formation in the western Mediterranean Sea: A high-resolution three-dimensional numerical study. *Dynamics of Atmospheres and Oceans*, 15(3–5), 301–332. https://doi.org/10.1016/0377-0265(91)90024-a

Madec, G., Delecluse, P., Crepon, M., & Chartier, M. (1991). A three-dimensional numerical study of deep-water formation in the northwestern Mediterranean Sea. *Journal of Physical Oceanography*, 21(9), 1349–1371. https://doi.org/10.1175/1520-0485(1991)021(1349:atdnso)2.0.co;2

Margirier, F., Testor, P., Heslop, E., Mallil, K., Bosse, A., Houpert, L., et al. (2020). Abrupt warming and salinification of intermediate waters interplays with decline of deep convection in the northwestern Mediterranean Sea. *Scientific Reports*, 10(1), 20923. https://doi.org/10.1038/s41598-020-77859-5

Marshall, J., & Schott, F. (1999). Open-ocean convection: Observations, theory, and models. *Reviews of Geophysics*, 37(1), 1–64. https://doi.org/10.1029/98rg02739

MEDOC. (1970). Observation of formation of deep water in the mediterranean sea, 1969. *Nature*, 227(5262), 1037–1040. https://doi.org/10.1038/2271037a0

Mertens, C., & Schott, F. (1998). Interannual variability of deep-water formation in the northwestern Mediterranean. *Journal of Physical Oceanography*, 28(7), 1410–1424. https://doi.org/10.1175/1520-0485(1998)028(1410:ivodwf)2.0.co;2

Millot, C., & Taupier-Letage, I. (2005). Circulation in the Mediterranean Sea. In *The Mediterranean Sea* (pp. 29–66). Springer Berlin Heidelberg. https://doi.org/10.1007/b107143

Napolitano, E., Iacono, R., Ciuffardi, T., Reseghetti, F., Poulain, P.-M., & Notarstefano, G. (2019). The Tyrrhenian intermediate water (TIW): Characterization and formation mechanisms. *Progress in Oceanography*, 170, 53–68. https://doi.org/10.1016/j.pocean.2018.10.017

NASA. (2023). QuikSCAT surface wind fields [Dataset]. Retrieved from https://cersat.ifremer.fr/Data/Latest-products/Reprocessing-Daily-QuikSCAT-Surface-Wind-Fields
Parras-Berrocal, I. M., Vázquez, R., Cabos, W., Sein, D. V., Álvarez, O., Bruno, M., & Izquierdo, A. (2022). Surface and intermediate water changes triggering the future collapse of deep water formation in the north western Mediterranean. *Geophysical Research Letters*, 49(4). https://doi.org/10.1029/2021gl095404

Quentin-Boris, B., Bruno, Z., Franck, D., Pierre, G., & Jonathan, B. (2020). Analysis of specific water masses transports in the western Mediterranean in the medrys1v2 twenty-one-year reanalysis [Poster]. FRANCE. Retrieved from https://archimer.ifremer.fr/doc/00628/74028/

Robinson, A., Leslie, W., Theocharis, A., & Lascaratos, A. (2001). Mediterranean Sea circulation. In *Encyclopedia of ocean sciences* (pp. 1689–1705). Elsevier. https://doi.org/10.1006/rwos.2001.0376

Ruti, P. M., Somot, S., Giorgi, F., Dubois, C., Flaounas, E., Obermann, A., et al. (2016). Med-CORDEX initiative for Mediterranean climate studies. *Bulletin of the American Meteorological Society*, 97(7), 1187–1208. https://doi.org/10.1175/bams-d-14-00176.1

Send, U., & Testor, P. (2017). Direct observations reveal the deep circulation of the western Mediterranean Sea. *Journal of Geophysical Research: Oceans*, 122(12), 10091–10098. https://doi.org/10.1002/2016jc012679

Severin, T., Kessouri, F., Rembauville, M., Sánchez-Pérez, E. D., Oriol, L., Caparros, J., et al. (2017). Open-ocean convection process: A driver of the winter nutrient supply and the spring phytoplankton distribution in the northwestern Mediterranean Sea. *Journal of Geophysical Research: Oceans*, 122(6), 4587–4601. https://doi.org/10.1002/2016jc012664

Skamarock, W., Klemp, J., Dudhia, J., Gill, D., Barker, D., Wang, W., et al. (2008). A description of the advanced research WRF version 3 (Tech. Rep.). https://doi.org/10.5065/D68S4MVH

Somot, S. (2005). *Modélisation Climatique du Bassin Méditerranéen: Variabilité et Scénarios de Changement Climatique (phdthesis)*. Université Toulouse III.

Somot, S., Houpert, L., Sevault, F., Testor, P., Bosse, A., Taupier-Letage, I., et al. (2016). Characterizing, modelling and understanding the climate variability of the deep water formation in the north-western Mediterranean Sea. *Climate Dynamics*, 51(3), 1179–1210. https://doi.org/10.1007/s00382-016-3295-0

Somot, S., Sevault, F., & Déqué, M. (2006). Transient climate change scenario simulation of the Mediterranean Sea for the twenty-first century using a high-resolution ocean circulation model. *Climate Dynamics*, 27(7–8), 851–879. https://doi.org/10.1007/s00382-006-0167-z

Testor, P., Bosse, A., Houpert, L., Margirier, F., Mortier, L., Legoff, H., et al. (2018). Multiscale observations of deep convection in the northwestern Mediterranean Sea during winter 2012–2013 using multiple platforms. *Journal of Geophysical Research: Oceans*, 123(3), 1745–1776. https://doi.org/10.1002/2016jc012671

Testor, P., & Gascard, J.-C. (2003). Large-scale spreading of deep waters in the western Mediterranean Sea by submesoscale coherent eddies. *Journal of Physical Oceanography*, 33(1), 75–87. https://doi.org/10.1175/1520-0485(2003)033(0075:lssodw)2.0.co;2

The Lab Sea Group. (1998). The Labrador Sea deep convection experiment. *Bulletin of the American Meteorological Society*, 79(10), 2033–2058. https://doi.org/10.1175/1520-0477(1998)079(2033:tlsdce)2.0.co;2

Turner, J. S. (1973). In U. of Cambridge (Ed.), *Buoyancy effects in fluids*. Cambridge University Press.

Waldman, R., Herrmann, M., Somot, S., Arsouze, T., Benshila, R., Bosse, A., et al. (2017a). Impact of the mesoscale dynamics on ocean deep convection: The 2012–2013 case study in the northwestern Mediterranean Sea. *Journal of Geophysical Research: Oceans*, 122(11), 8813–8840. https://doi.org/10.1002/2016jc012587

Waldman, R., Somot, S., Herrmann, M., Bosse, A., Caniaux, G., Estournel, C., et al. (2017b). Modeling the intense 2012–2013 dense water formation event in the northwestern Mediterranean Sea: Evaluation with an ensemble simulation approach. *Journal of Geophysical Research: Oceans*, 122(2), 1297–1324. https://doi.org/10.1002/2016jc012437

Waldman, R., Somot, S., Herrmann, M., Sevault, F., & Isachsen, P. E. (2018). On the chaotic variability of deep convection in the Mediterranean Sea. *Geophysical Research Letters*, 45(5), 2433–2443. https://doi.org/10.1002/2017gl076319