# Near-bed stratification controls bottom hypoxia in ice-covered alpine lakes

Marie-elodie Perga, Camille Minaudo, Tomy Doda, Florent Arthaud, Harsh Beria, Hannah E Chmiel, Nicolas Escoffier, Thibault Lambert, Raphaelle Napolleoni, Biel Obrador, et al.

► To cite this version:

Marie-elodie Perga, Camille Minaudo, Tomy Doda, Florent Arthaud, Harsh Beria, et al. Near-bed stratification controls bottom hypoxia in ice-covered alpine lakes. Limnology and Oceanography, 2023, 68 (6), pp.1232-1246. 10.1002/lno.12341. hal-04222315

HAL Id: hal-04222315

https://hal.inrae.fr/hal-04222315

Submitted on 29 Sep 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons Attribution 4.0 International License
LIMNOLOGY
and
OCEANOGRAPHY

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

Linnol. Oceanogr. 68, 2023, 1232–1246
© 2023 The Authors. *Limnology and Oceanography* published by Wiley Periodicals LLC on
behalf of Association for the Sciences of Limnology and Oceanography.
doi: 10.1002/lno.12341

Near-bed stratification controls bottom hypoxia in ice-covered
alpine lakes

Marie-Elodie Perga ID, ¹,²* Camille Minaudo ID, ³,⁴ Tomy Doda ID, ⁵ Florent Arthaud, ⁶ Harsh Beria, ⁷

Hannah E. Chmiel,⁴ Nicolas Escoffier ID, ¹ Thibault Lambert, ¹ Raphaelle Napolleoni, ⁸ Biel Obrador ID, ⁹

Pascal Perolo ID, ¹ Janine Rüegg, ¹ Hugo Ulloa ID, ¹⁰ Damien Bouffard ID, ¹,⁵

¹Faculty of Geoscience and Environment, University of Lausanne, Lausanne, Switzerland

²Geopolis, UNIL (on behalf of the Lacs Sentinelles Network), Lausanne, Switzerland

³Department of Evolutionary Biology, Ecology and Environmental Sciences, University of Barcelona, Barcelona, Spain

⁴Physics of Aquatic Systems Laboratory, Margaretha Kamprad Chair, Institute of Environmental Engineering, École Polytechnique Fédérale de Lausanne, Lausanne, Switzerland

⁵Department of Surface Waters—Research and Management, Eawag, Swiss Federal Institute of Aquatic Science and Technology, Kastanienbaum, Switzerland

⁶UMR CARRTEL, University Savoie-MontBlanc, Le Bourget du Lac, France

⁷Department of Environmental Systems Science, ETH Zurich, Zurich, Switzerland

⁸ASTERS (on behalf of the Lacs Sentinelles Network), Pringy, France

⁹Institut de Recerca de la Biodiversitat (IRBio), University of Barcelona, Barcelona, Spain

¹⁰Department of Earth and Environmental Science, University of Pennsylvania, Philadelphia, Pennsylvania

Abstract

In ice-covered lakes, near-bottom oxygen concentration decreases for most of the wintertime, sometimes down to the point that bottom waters become hypoxic. Studies insofar have reached divergent conclusions on whether climate change limits or reinforces the extent and duration of hypoxia under ice, raising the need for a comprehensive understanding of the drivers of the dissolved oxygen (DO) dynamics under lake ice. Using high-temporal resolution time series of DO concentration and temperature across 14 mountain lakes, we showed that the duration of bottom hypoxia under ice varies from 0 to 236 d within lakes and among years. The variability of hypoxia duration was primarily explained by changes in the decay rate of DO above the lake bottom rather than by differences in DO concentration at the ice onset or in the ice-cover duration. We observed that the DO decay rate was primarily linked to physical controls (i.e., deep-water warming) rather than biogeochemical drivers (i.e., proxies for lake or catchment productivity). Using a simple numerical model, we provided a proof-of-concept that the near-bed stratification can be the mechanism tying the DO decay rate to the sediment heat release under the ice. We ultimately showed that the DO decay rate and hypoxia duration are driven by the summer light climate, with faster oxygen decline found under the ice of clearer cryostratified alpine lakes. We derived a framework theorizing how the hypoxia duration might change under the ice of alpine lakes in a warmer climate.

This is an open access article under the terms of the [Creative Commons Attribution-NonCommercial License](https://creativecommons.org/licenses/by-nc), which permits use, distribution and reproduction in any medium, provided the original work is properly cited and is not used for commercial purposes.

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

**Author Contribution Statement:** Marie-Elodie Perga and Camille Minaudo designed the study, all authors contributed to the data analysis during the Lakathon of February 2020. Damien Bouffard, Hugo N. Ulloa, and Tomy Doda designed the testbed model. Marie-Elodie Perga and Camille Minaudo drafted the manuscript, all authors contributed to the final version of the manuscript.

The cryosphere is shrinking at well-documented rates, but the far-reaching ecosystemic implications of ice loss remain poorly anticipated (Fountain et al. 2012). This holds for cryosphere features such as sea ice, snow, and glaciers, with lake ice making no exception. Current changes in the phenology and extent of lake ice have been regularly updated and refined since the seminal work by Magnuson et al. (2000). A recent field-based study estimated that the ice cover of Northern hemisphere lakes has, on average, shortened by 17 d over the last century (Sharma et al. 2021). Significant achievements in remote sensing allow the monitoring of lake ice phenology where ground data are scarce (Wang et al. 2021). Deterministic or data-driven modeling approaches can predict the vulnerability of currently dimictic lakes to permanent ice loss until the end of

<footnote>*Correspondence: marie-elodie.perga@unil.ch</footnote>
the century (Woolway and Merchant 2019; Caldwell et al. 2021) as well as more continuous changes in ice phenology (Råman Vinnå et al. 2021). While the modifications of winter ice conditions and phenology are expected to result in major rearrangements of physical and ecological processes in lakes losing ice cover, most still have to be evidenced and quantified (Salonen et al. 2009; Hampton et al. 2017; Cavaliere et al. 2021).

Pioneering investigations of the processes occurring under lake ice had to build on a limited number of observations due to the difficulty of sampling (Mortimer and Mackereth 1958). The development of sensor technology combined with a renewed interest in wintertime limnology has generated large amounts of high-resolution data over the last decade. However, rather than clearing up the picture, recent studies have unveiled diverse and complex physical, biogeochemical, and ecological processes occurring under the ice (Obertegger et al. 2017; Song et al. 2019; Zdorovennova et al. 2021). Conceptual frameworks have been developed to unify pioneering and more recent knowledge of under-ice processes and derive assumptions on the ecological consequences of a changing ice cover on lakes (Cavaliere et al. 2021; Jansen et al. 2021). All conceptual frameworks place dissolved oxygen (DO) concentrations under the ice as a central variable both integrating the ecological responses to changing wintertime conditions and determining the effects of ice phenology alteration on the following seasons (Cavaliere et al. 2021; Jansen et al. 2021).

Processes acting upon DO under the ice can include a significant amount of primary production, depending on the thickness and optical properties of the ice (Zdorovennova et al. 2021). The DO stock under the ice nevertheless diminishes with time due to the lake surface isolation from the atmosphere and the predominance of heterotrophy (Leppi et al. 2016; Obertegger et al. 2017; Song et al. 2019). Therefore, in many ice-covered lakes, the near-bottom oxygen concentration decreases for most of the wintertime, sometimes down to the point that bottom waters become hypoxic (Mathias and Barica 1980; Malm et al. 1998; Stefanovic and Stefan 2002). The vertical extent of the hypoxic layer and duration of hypoxia determine the availability of oxic refuges for many organisms actively overwintering under the ice (Davis et al. 2020) and, over the long-term, shape biological communities according to their abilities and adaptations to low ambient DO conditions (Magnuson et al. 1985). With oxygen depletion, microbial metabolism shifts to other electron-acceptors, modifying the partitioning and concentrations of major elements (N, C, S, P, Fe) under the ice, with potentially far-reaching consequences for the following open water season (Gammons et al. 2014; Kincaid et al. 2022). While climate warming is expected to result in significant modifications of the DO dynamics under lake ice, studies insofar have reached divergent conclusions on whether climate change shall limit or instead reinforce the duration and vertical extent of hypoxia under the ice. For instance,

climate warming was observed or expected to promote oxygenation under ice in some cases, (Flaim et al. 2020; Smits et al. 2021) but predicted to extend the duration of hypoxia under the ice for others (Golosov et al. 2007; Couture et al. 2015).

Schematically, the dynamics and duration of bottom hypoxia under lake ice depend on three variables: the ice-cover duration, the DO stock at the ice onset, and the rate of DO decay under the ice (Terzhevik et al. 2009); all of which are susceptible and expected to vary with climate warming. A shorter ice cover reduces the time over which hypoxia can develop and persists (Smits et al. 2021). Besides, climate warming affects water temperatures and mixing duration in fall (Råman Vinnå et al. 2021), both influencing DO stocks at the ice onset (Flaim et al. 2020). The DO decay rates are also likely to respond to climate change, potentially through biogeochemical but also physical pathways. For instance, DO decay under the ice primarily arises from the sediment oxygen uptake (SOU; Bengtsson and Ali-Maher 2020), i.e., the diffusion of oxygen within the sediment where it is consumed by both the aerobic mineralization of sediment organic matter (OM) and the abiotic oxidation of reduced solutes produced by anaerobic microbial metabolism (Müller et al. 2012). The oxygen consumption in the sediment increases with the sediment content in OM (Müller et al. 2012), so that a greater SOU is expected in situations where climate warming results in greater OM exports to the sediment (Steinsberger et al. 2017), either through enhanced lake primary production (Preston et al. 2016) or greater terrestrial OM inputs (Couture et al. 2015). Yet, DO decay might not be solely biogeochemically controlled. In the first investigations of under-ice dynamics, DO decay rates were already observed to increase with the warming rate of deep waters by sediment heat release (Mortimer and Mackereth 1958), an association confirmed by more recent studies (Golosov et al. 2007; Terzhevik et al. 2009). Although the processes tying the sediment heat release to DO decay remain to be identified, sediment heat content is another pathway by which climate warming might modify DO decay under lake ice (Fang and Stefan 1996; Golosov et al. 2007).

Because of the diversity of involved processes contributing to the DO dynamics and the climate factors acting upon them, there is no single scenario for the effect of climate warming on DO under lake ice. The different components of the DO dynamics can respond synchronously and sometimes through divergent pathways to climate warming. For instance, a modeling approach applied to the boreal Lake Langtjern suggested that despite an earlier ice break up, DO under ice is likely to decrease due to greater climate-driven dissolved organic carbon inputs from the catchment (Couture et al. 2015). Therefore, shifts in the extent and degree of hypoxia in ice-covered lakes due to climate change are case-specific.

Processes of DO decay in high altitude lakes, where deep-water hypoxia can develop over the course of the 6–9 months
of ice cover despite low productivity, are even less understood (Smits et al. 2021). Compared to the low elevation temperate and boreal lakes, their high-altitude counterparts are less prone to summer hypoxia since the absence or short duration of summer stratification limits deep deoxygenation during the open water season (Perga et al. 2018; Granados et al. 2020). However, high-altitude lakes are exposed to amplified climate warming (Gobiet et al. 2014) to which they commonly respond by a later ice-on and an earlier ice-off (Parker et al. 2008; Brown and Duguay 2010; Preston et al. 2016; Smits et al. 2021), and possibly by a longer fall mixing (Flaim et al. 2020). The consequences of climate warming on hypoxia duration during wintertime remain unclear as the outcome would depend on the relative importance of ice duration, DO stocks at ice onset and DO decay rates in controlling near-bottom DO under the ice.

Here we present a comprehensive database of high-
temporal resolution time series of DO concentration and tem-
perature across 14 mountain lakes located in the French Alps,
joined with a set of attributes for the lakes and their catch-
ment. First, we quantified the duration of hypoxia spatially,
i.e., among lakes, and temporally, i.e., among years within
lakes. We then identified which of ice phenology, DO concen-
trations at the ice onset or DO decay rate was the main driver
of hypoxia duration. As DO decay rate stood out as the domi-
nant predictor of hypoxia duration, we evaluated, in a second
step, the relative contribution of the biogeochemical and
physical controls to DO decay rates. In the final step, we pro-
posed a simple mechanistic model as a proof-of-concept of the
mechanisms derived from the data analysis.

**Materials and methods**

Study sites and monitoring setup

The studied lakes are part of the observatory of the Sentinel
lakes network (http://www.lacs-sentinelles.org/, last accessed
February 2021, Birck et al. 2013). The network developed a set
of standardized protocols adapted for long-term monitoring of
mountain lakes, accounting for their remoteness and extreme
conditions. The protocols rely on high-frequency monitoring
from autonomous sensors and one annual field campaign in
early September to maintain the sensors, download the data,
and perform vertical profiling, Secchi depth measurements,
and water sampling. Fourteen of these lakes are considered in
the analysis and are classified as dimictic or cold monomictic
(Fig. 1, Table 1).

The monitored lakes span the French Alps, from the foot-
hill of the Mont Blanc down to the Southernmost part of the
Alps (Fig. 1), with altitudes ranging from 1945 to 2752 m
above sea level (Table 1). They are small lakes with surface
areas and maximum depths (zmax) that vary over one order of
magnitude (from 1 to 16 ha, 7 to 50 m deep). Drainage ratios
(catchment to lake surface areas) range from 9 to 130. All
lakes are oligotrophic, poorly productive, and ice-covered for

150–230 d of the year. Despite a narrow range of altitudes, the
lake catchments are highly variable in terms of vegetation cov-
erage (Normalized Difference Vegetation Index [NDVI] from
0 to 0.7) and proportions of cryogenic features (glaciers and
permafrost from 0% to 73% and 0% to 94%, respectively).

As part of the monitoring protocols, all lakes are equipped
with autonomous temperature sensors (Tinytag® or Hobo
Water Temperature Pro v2®, accuracy ± 0.5 °C) 2 m below the
surface, along with a temperature-DO sensor (PME Minidot,
accuracy ± 0.1°C and ± 0.3 mg O₂ L⁻¹) placed 1 m above the
sediment and anchored to a submerged buoy at the maximum
lake depth. Both water temperature and DO concentration
were recorded at 1-h intervals. Twenty-one DO and water tem-
perature data series were available.

DO dynamics

Since no direct ice measurements were available, the ice-
cover phenology was estimated from high-frequency near-
surface and near-bottom water temperatures. The date of
ice-on was defined as the time when lake bottom tempera-
ture reached a local minimum after the fall cooling,
i.e., the moment when the water column suddenly isolates
from the atmosphere and the inverse thermal stratifica-
tion, resulting from a coupled effect of heat loss and
sediment heat release, becomes a distinctive feature of
under-ice waters (Fang and Stefan 1996). The date of ice-
off was defined as the time when the near-surface water
temperature was rising after ice-on by more than 0.1°C
day⁻¹, i.e., the moment when the water column gets in
contact with the atmosphere and warms up (SI 1). We
acknowledge that the quantification of the ice cover phe-
nology is subject to arbitrary choices but keeping consis-
tent metrics throughout the study ensures that trends are
conservative among sites and years.

DO\_init was computed as the near-bottom DO at the deter-
mined date of ice on. Time to hypoxia corresponded to the
number of days to reach DO = 2 mg L⁻¹ after the date of the
ice onset. The hypoxia duration was computed as the number
of days for which the near-bottom DO ≤ 2 mg L⁻¹ after the ice
onset and before the determined ice-off.

In earlier studies, under-ice DO depletion had been
modeled from either a zero- or a first-order kinetics. Zero-order
kinetics model the DO dynamics as a linear decrease over
time, computed over the first weeks after ice onset or in-
between changepoints in the DO time-series (Obertegger
et al. 2017; Granados et al. 2020; Smits et al. 2021). A zero-
order model assumes a substrate limitation of DO consump-
tion which, close to the lake bottom (Beutel 2003), might be
better adapted to model the oxidation of reduced compounds
released from the sediment (Müller et al. 2012). In contrast,
first-order kinetics (Golosov et al., 2007, Terzhevik et al. 2009)
model the DO dynamics as an exponential decay that better
reflects diffusion-limited processes at the sediment–water
interface (Beutel 2003), i.e.,
<figure><img src="image_3.png" /><figcaption>Fig. 1. Locations of the 14 monitored alpine lakes, distributed across the French Alps.</figcaption></figure>

$$DO(t) = DO_{init} \times e^{-K_{exp}}(t-t_{0})$$

where $t_0$ is the determined ice-on ($t$ is expressed in days), $DO_{init}$ is the initial DO concentration. $K_{exp}$ is the exponential decay coefficient that is fitted to the data. We opted for the exponential decay model as it fitted better, was the most consistent across lakes, but also over time as it reproduces the asymptotic trend as DO gets closer to 0 mg L⁻¹. Conclusions were however independent of the type of model used to simulate the DO dynamics under ice.

In several cases, near-bottom DO first increased for several days after ice-on, likely due to residual light penetration through the ice allowing either under-ice primary production (Granados et al. 2020), or due to buoyancy-driven downslope flows (Ulloa et al. 2019). To minimize this bias, the decay models were fitted on the DO time series starting after the local DO maximum was reached. The length of the fitting window was manually adapted

in the case of two lakes, where the observed DO patterns showed complex features. This was the case in Lake Bresses Inférieur during the winter of 2018–2019, where DO followed an exponential decay for 50 d after ice-on with significant unexplained upward deviations afterwards. The curve-fitting was conducted over the first 50 d after ice-on in this case. In Lake Merlet Supérieur, DO measurements presented an abrupt increase 100 d after ice-on and $K_{exp}$ was fitted to the data recorded before this event. Overall, fitting performances were statistically robust, with determination coefficients $r^2>0.83$, and a mean $r^2$ of 0.93 across all lakes and years.

## Explanatory variables

We investigated the relationships existing between DO dynamics under ice and potential predictors describing lake morphology, thermal stratification at ice-onset (i.e., cryoregime), lake productivity, and the warming rate of deep
<table><caption>**Table 1.** Main geographical attributes for 14 alpine lakes and their respective catchments.</caption><thead><tr><th>Lake</th><th>Mixing regime</th><th>Altitude [m a.s.l.]</th><th>Area [ha]</th><th>z<sub>max</sub> [m]</th><th>Secchi depth [m]</th><th>Drainage ratio [unitless]</th><th>Catchment NDVI [dl]</th><th>Macrophyte coverage [%]</th><th>Monitored winters (n)</th></tr></thead><tbody><tr><td>Anterne</td><td>Dimictic</td><td>2063</td><td>11.1</td><td>12</td><td>3.6–10.0</td><td>21</td><td>0.60</td><td>9</td><td>2015–2019 (4)</td></tr><tr><td>Arpont</td><td>Cold monomictic</td><td>2666</td><td>3.7</td><td>17</td><td>0.5–0.9</td><td>91</td><td>0.02</td><td>ND</td><td>2018–2019 (1)</td></tr><tr><td>Bresses inférieur</td><td>Dimictic</td><td>2458</td><td>1.2</td><td>13</td><td>6.0–9.0</td><td>9</td><td>0.35</td><td>ND</td><td>2018–2019 (1)</td></tr><tr><td>Bresses supérieur</td><td>Dimictic</td><td>2501</td><td>1</td><td>12</td><td>5.0–7.0</td><td>36</td><td>0.20</td><td>ND</td><td>2018–2019 (1)</td></tr><tr><td>Cornu</td><td>Dimictic</td><td>2276</td><td>6.4</td><td>22</td><td>10.0–13.0</td><td>10</td><td>0.37</td><td>0</td><td>2018–2019 (1)</td></tr><tr><td>Lauzanier</td><td>Cold monomictic</td><td>2284</td><td>3.8</td><td>7</td><td>7.0</td><td>130</td><td>0.34</td><td>30</td><td>2017–2019 (2)</td></tr><tr><td>Merlet</td><td>Dimictic</td><td>2452</td><td>4.9</td><td>29</td><td>7.6–13.3</td><td>13</td><td>0.40</td><td>6</td><td>2018–2019 (1)</td></tr><tr><td>Mont Coua</td><td>Cold monomictic</td><td>2672</td><td>2.5</td><td>10</td><td>3.8–4.1</td><td>32</td><td>0.14</td><td>ND</td><td>2018–2019 (1)</td></tr><tr><td>Muzelle</td><td>Dimictic</td><td>2105</td><td>9.7</td><td>18</td><td>0.6–8.0</td><td>46</td><td>0.24</td><td>&lt;1</td><td>2014–2015;<br>2018–2019</td></tr><tr><td>Noir du Carro</td><td>Cold monomictic</td><td>2759</td><td>1.2</td><td>11</td><td>10.0–11.3</td><td>41</td><td>0.16</td><td>0</td><td>2018–2019 (1)</td></tr><tr><td>Petarel</td><td>Dimictic</td><td>2095</td><td>2.5</td><td>20</td><td>11.6–13.5</td><td>53</td><td>0.44</td><td>0</td><td>2018–2019 (1)</td></tr><tr><td>Pisses</td><td>Cold monomictic</td><td>2490</td><td>1.7</td><td>7</td><td>7.5–8.0</td><td>47</td><td>0.44</td><td>ND</td><td>2017–2019 (2)</td></tr><tr><td>Pormenaz</td><td>Dimictic</td><td>1945</td><td>4.4</td><td>9</td><td>5.7–8.4</td><td>11</td><td>0.71</td><td>10</td><td>2018–2019 (1)</td></tr><tr><td>Rabuons</td><td>Dimictic</td><td>2500</td><td>16.3</td><td>50</td><td>11.5–13.5</td><td>19</td><td>0.19</td><td>&lt;1</td><td>2017–2019 (2)</td></tr></tbody></table>

waters. Morphological descriptors were lake depth and lake photic ratio; the latter being defined as the ratio of light penetration depth in water to the maximum lake depth. The lake photic ratio (*phot ratio*) was compiled for each year such as,

$$phot_{ratio} = \frac{z_{eu}}{z_{max}}$$

with the euphotic depth $z_{eu} = 2.5 \cdot z_{Secchi}$ (Joint Research Centre, Institute for Environment Sustainability 2012).

We characterized the lake cryoregime (i.e., cryomictic or cryostratified) as defined by Yang et al. (2021) from the mean temperature of the water column at ice onset (*WTR ice-on*). Briefly, cryomictic lakes are shallow or wind-exposed lakes in which waters can mix and thus cool down to 0–1°C before ice forms. In cryostratified lakes, ice forms as soon as the surface reaches the freezing point, while deeper waters are 2–4°C (Yang et al. 2021). Lakes with *WTR ice-on* in the range 0–1°C were classified as cryomictic, whereas lakes with *WTR ice-on* above 2°C were classified as cryostratified (Yang et al. 2021).

Biogeochemically, the DO decay rate is expected to be greater in lakes with greater supplies of organic carbon (Leppi et al. 2016; Terry et al. 2017) from both autochthonous and terrigenous origins. In alpine lakes, Secchi depth alone is a poor predictor of pelagic productivity because water transparency strongly varies with the catchment sediment load (Olson et al. 2017; Perga et al. 2018). Considering that primary production is essentially littoral in small and clearwater oligotrophic lakes (Devlin et al. 2016), the percent macrophytes coverage was considered as a proxy for lake productivity. It was measured in August 2016 or 2017 in 25 quadrats regularly distributed on the lake area. Inputs of terrigenous OM and

catchment productivity were estimated using drainage ratios (Cremona et al. 2019), and the catchment vegetation cover (NDVI; Larsen et al. 2011), respectively. The NDVI values correspond to the mean peak productivity over 2013–2018.

The warming rate of deep waters was assessed through its e-folding timescale coefficient ($K_{heat}$), computed by fitting the near-bottom daily-averaged temperature measurements during the ice-covered period to the model $T_t = T_f - (T_f - T_0)e^{-K_{heat}(t-t_0)}$, where $T_0$ corresponds to the temperature measured at the time of ice-on, $T_f$ at ice melt, and $t$ is the time expressed in days after ice-on ($t_0$).

The relationships between DO decay rate coefficient, i.e., $K_{exp}$, and potential predictors were tested using linear regression models, with log-normalized data when necessary. Some lakes had records for multiple winters (Anterne, Lauzanier, Muzelle, Pisses, Rabuons), while some predictors (lake morphology, catchment land cover) were year-invariant. Weighting the contribution of lakes by the number of years of measurements ensured a balanced contribution of lakes to the final models. However, results did not vary between weighted and unweighted models.

### Reaction-diffusion model

To assess the mechanistic processes linking predictors to DO decay in the field dataset, we used a simplified model of coupled DO and water temperature under the ice, set on a hypothetical alpine lake with a parabolic morphology of 20 m max depth (average max depth was 17 m in our dataset) and 5 ha surface area (average surface area was 4.8 m in our dataset). The modeling framework considers diffusion-reaction equations for temperature and oxygen. The model
solves the one-dimensional (vertical), diffusion reaction equation for horizontally averaged scalar field, C(z,t) being either the DO or the T field

$$A(z) \frac{\partial C(z, t)}{\partial t} = \frac{\partial}{\partial z} \left( A(z) K(z, t) \frac{\partial C(z, t)}{\partial z} \right) + \frac{\partial A(z)}{\partial z} F$$

with,

* A(z) the bathymetry;

* K(z,t) the vertical turbulent diffusivity defined as a function of the buoyancy frequency (N²): Kz = aN⁻²b (Jassby and Powell 1975; Heinz et al. 1990; Hondzo and Stefan 1993);

* F is the flux at slopping bottom (∂A/∂z). F is either the sediment heat flux or the sediment oxygen uptake:

* Case DO: F = -kDO (DO(z) - DOswi), with kDO = DDO/δ, and DOswi = DDO(z)/2 with DDO the molecular diffusivity for DO, δ the constant diffusive boundary layer and DOswi the oxygen concentration at the sediment water interface (Müller et al. 2012; Bouffard et al. 2013).

* Case T: F = -kT, with kT = He^kTz e^(-Tt), with ρ, the water density, Cp, the heat capacity for water, τ, the decay time scale for sediment heat flux, kdl, the light extinction and H the sediment heat flux per unit area. Here the sediment heat flux decreases over time and depth.

The ice-water boundary is treated as a Dirichlet condition for temperature, T(t,z=0) = 0°C, and as a no flux Neuman condition for DO, (dC/dz)|z=0 = 0). The lake bottom (i.e., deepest point) boundary condition is treated as no-flux Neuman condition for both T and DO (the flux at the bottom is already included in the last term of the equation). We stress that the equation is valid in the case of horizontally averaged constituents (Stepanenko et al. 2016). Full model specifications and parameters are provided in SI 2.

At ice onset, the DO concentration was considered uniform over depth, whereas the initial vertical temperature profile was determined by the respective cryoregime. The near-bottom (1 m above the sediment) DO decay rate was computed for different scenarios, testing for the effects of:

i. Cryoregimes, i.e., a cryomictic (well mixed water column at ice onset with temperature close to 0°C with a temperature increasing near-bed) and a cryostratified regime (temperature increase with depth close to the surface and well mixed conditions down to the bottom; Yang et al. 2021). Boundary conditions for water temperature were set at 0°C at the surface and 3°C at the lake bottom.

ii. Initial sediment heat release set at 0, 1, 2, and 4 W m⁻². Values for the initial sediment heat flux were determined based on simulations on a well-established bulk model FLake (Golosov and Kirillin 2010) for a lake of similar latitude and longitude as Lake Anterne with a light extinction coefficient ranging between 0.4 and 4 m⁻¹ (i.e., Secchi depth within 0.5–5 m, comparable to the range found in

our dataset). Such values for the initial sediment heat release are in line with those reported in previous studies (Kirillin et al. 2012; Zdorovennova et al. 2021)

## Results

### Variability of under-ice DO dynamics across and within lakes

In all lakes, bottom DO decreased under ice, with total DO losses for the ice-on period ranging from 2.2 to 12.3 mg L⁻¹ for Lakes Arpont in winter 2018–2019 and Pisses in 2017–2018, respectively. The duration of hypoxia was highly variable among lakes, from no hypoxia for Lakes Arpont 2018–2019, Muzelle 2018–2019 and Rabuons 2018–2019, to 236 hypoxic days for Mont Coua, where hypoxia occurred as early as 24 d after the ice onset (Figs. 2, 3a) and persisted a month beyond the ice-off. The duration of hypoxia was also highly variable among years for a single lake. For instance, Lake Muzelle experienced 0 d of hypoxia in 2018–2019 and 108 d in 2014–2015. In Lake Anterne, hypoxia duration ranged between 143 and 195 d over the 4 recorded winters.

The duration of hypoxia was neither correlated to ice cover duration (Fig. 3a), nor to DO saturation at the ice onset (DOinit, Fig. 3b). However, the duration of hypoxia among lakes and years was strongly dependent on the DO decay rate Kexp (Fig. 3c). Because the coefficient Kexp was, by far, the best predictor of hypoxia duration in our dataset, the subsequent analyses focused on uncovering the processes governing the DO decay rate in mountain lakes.

The DO decay coefficient Kexp ranged from 0.0016 to 0.0626 d⁻¹, i.e., e-folding decay timescales between 598 and 16 d, respectively. Kexp showed substantial variability among lakes, and among years for a single lake (Fig. 4a). Some lakes showed temporary breaks in the DO decay trends (for instance, Merlet supérieur, Muzelle in 2014–2015 and Bresses inférieur, Fig. 2). In the case of Bresses inférieur, hypoxia was reached 100 d later than expected from an exponential decay model (Fig. 4b). For all other lakes, the DO dynamics fell within the expected model describing time to hypoxia from Kexp, confirming that an exponential decay model explains the temporal evolution of DO in most lakes over the entire ice-covered period (Fig. 4b).

All lakes, except for Lakes Mont Coua and Lauzanier 2018–2019, had bottom water temperatures above 2°C at ice onset (Fig. 5), indicating that the majority of lakes in our dataset can be regarded as cryostratified (Yang et al. 2021). Lake Mont Coua, with bottom water temperature as low as 0.2°C at the ice onset, can unambiguously be classified as cryomictic. For Lauzanier, the cryoregime is somewhat uncertain as bottom water temperature were slightly lower than 2°C in 2018–2019, but not as low as could be expected for a strictly cryomictic lake (Yang et al. 2021). Nevertheless, the greatest Kexp was observed in the single
unambiguously cryomictic lake (Mont Coua, Fig. 5). With
only one cryomictic lake in the dataset, we cannot robustly
conclude on the potential effect of cryoregime on DO decay
*K*<sub>exp</sub>. As the cryomictic lake appeared as an outlier in our
dataset, we analyzed the predictive power of environmental
variables for cryostratified lakes only. Note that when the
dataset was restricted to cryostratified lakes only (including
Lauzanier 2018–2019), *K*<sub>exp</sub> remained the best predictor for
hypoxia duration (*r*<sup>2</sup> = 0.67, *p*<10<sup>-5</sup>), far above initial DO
(*r*<sup>2</sup> = 0.11, *p* = 0.24) and ice duration (*r*<sup>2</sup> = 0.08, *p* = 0.23).

Predictors of under-ice DO decay rates for cryostratified lakes

Figure 6 summarizes the predictors of under-ice DO decay rates for cryostratified lakes. *K*<sub>exp</sub> was negatively correlated with lake depth (Fig. 6a), although the determination coefficient from linear regression was relatively low (28%, Table 2). *K*<sub>exp</sub> increased significantly with the photic ratio at the end of summer (Fig. 6b), following a linear relationship and providing the greatest coefficient of determination (59%, Table 2). *K*<sub>exp</sub> was positively correlated to the deep warming rate under the ice (*K*<sub>heat</sub>), which was a better predictor of *K*<sub>exp</sub> (Fig. 6c) than any of the variables used as OM supply proxies (i.e., drainage ratio, NDVI and macrophyte

cover, Fig. 6d–f). *K*<sub>heat</sub> explained up to 33% of *K*<sub>exp</sub>
variability.

Diffusion–reaction model simulations

Field data suggested that the lake cryoregime and the
e-folding coefficient for deep warming, i.e., *K*<sub>heat</sub>, were likely
to affect *K*<sub>exp</sub>, which we investigated further using model sim-
ulations. Specifically, the model considers the initial
cryoregime and resolves the vertical temperature and DO evo-
lution based on heat release from sediment interface and irre-
versible vertical diffusive fluxes across the water column.
Reported *K*<sub>exp</sub> are those simulated 1 m above the sediment. In
this hypothetical lake mechanistic model, the cryoregime had
a significant influence on *K*<sub>exp</sub>, with values about 25% higher
in a cryomictic lake compared to a cryostratified lake, in the
absence of any sediment heat release (Fig. 7a). *K*<sub>exp</sub> increased
along with a greater sediment heat release in the cryostratified
regime, while *K*<sub>exp</sub> sensitivity to sediment heat flux was
extremely low in the cryomictic lake (Fig. 7a). The vertical
temperature profile revealed that the near-bed stratification,
that already existed at the ice onset, persisted under the ice of
the cryomictic lake, even in the absence of sediment heat
release (Fig. 7b). Sediment uptake drew down DO from the
bottom layer only (up to 2–3 m above the sediment), and, as

<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Patterns of near-bottom (1 m above the sediment) DO concentration, shown both as measured (thin black lines) and fitted (thick red lines), during the under-ice period in the 14 studied lakes. The darker part of the thin lines indicates the data used for curve-fitting and computation of *K*<sub>exp</sub> coefficients. The threshold for hypoxia (2 mg L<sup>-1</sup>) is indicated by the horizontal gray line.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Relationships of bottom hypoxia duration under the ice and (a) the ice cover duration, (b) the initial near-bottom DO at ice-onset, and (c) the DO decay coefficient *K*<sub>exp</sub>. For Lake Mont Coua, the spring mixing occurred a month later than the ice melt, so that the hypoxia duration exceeds the ice cover duration (SI 1). The size of data points is proportional to the weight given within the linear model (see *Materials and methods*).</figcaption></figure>

the near-bed density gradient prevented DO diffusion and resupply from the above layers, resulted to a greater near-bottom DO depletion in a cryomictic regime (Fig. 7d). Because of the persistence of near-bed stratification, DO depletion was similar irrespective of the sediment heat flux for a cryomictic lake (Fig. 7a). In the cryostratified regime, the vertical temperature

profile remained homogenous below the 14 m depth in the absence of sediment heat release (Fig. 7c) so that the downward DO diffusion in the water column partly compensated for sediment uptake. The near-bottom DO decay was limited in the cryostratified regime when sediment heat release was absent (Fig. 7e). With greater sediment heat release, a near-bed

<figure><img src="image_6.png" /><figcaption>**Fig. 4.** (a) DO decay rate coefficient *K*<sub>exp</sub> for all lakes (x axis) and years (colors). (b) Relationship between observed *K*<sub>exp</sub> and time (in days) for lakes to reach near bottom hypoxic conditions, i.e., DO ≤ 2 mg L<sup>-1</sup>. Data points represent the observed time to hypoxia after the onset of ice in days. The dotted line represents the theoretical first-order differential equation model estimating the time to hypoxia τ<sub>hypoxia</sub> based on observed *K*<sub>exp</sub>, such that τ<sub>hypoxia</sub> = *K*<sub>exp</sub><sup>-1</sup> log(<sup>DO<sub>init</sub></sup>/<sub>DO<sub>hypox</sub></sub>), with DO<sub>hypox</sub> the threshold for hypoxia, set at 2 mg L<sup>-1</sup>, and assuming DO<sub>init</sub> = 9 mg L<sup>-1</sup> (the mean value of the dataset). For lakes that never reached hypoxia, the time to hypoxia was set at 230 d (the longest ice duration).</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Fig. 5.** DO decay rate $K_{exp}$ and cryoregimes classification based on mean water temperature in the water column at ice-on. Lake Mont Coua was the only unambiguously cryomictic lake in our dataset.</figcaption></figure>

temperature gradient developed under the ice of the cryostratified lake (Fig. 7c), reducing DO downward diffusion from the upper layers and leading to faster DO decay near the bottom (Fig. 7e).

## Discussion

Our extensive dataset revealed highly variable near-bottom DO dynamics under the ice of alpine lakes. The rate of DO decay $K_{exp}$ was the most critical parameter for development of near-bed hypoxia. The value of $K_{exp}$ allowed confident prediction of the time required to reach near-bottom hypoxia as well as of the duration of hypoxia under ice. The DO decay rate was twice as high for the cryomictic lake as compared to the other, cryostratified, lakes, suggesting that the lake cryoregime could influence the DO dynamics under the ice. The DO decay rate among cryostratified lakes varied mainly with light climate (i.e., the photic ratio). The statistical analyses

<figure><img src="image_8.png" /><figcaption>**Fig. 6.** Exponential decay coefficients $K_{exp}$ (per day) as a function of the potential predictors, i.e., (a) mean lake depth, (b) the photic ratio, (c) deep warming rate at the lake near-bottom, (d) lake drainage ratio, (e) catchment vegetation cover, NDVI, and (f) lake relative macrophyte coverage. Linear fits are performed only over cryostratified lakes (closed points), while the cryomictic lake (open point) data is indicated only for a comparative purpose. The size of data points is proportional to the weight given within the linear model.</figcaption></figure>
<table><caption>**Table 2.** Linear model results (determination coefficient r² and associated p-value) for K<sub>exp</sub> and the selected variables for the cryostratified lakes.</caption><thead><tr><th>Predictors</th><th>Variable (transformation; unit)</th><th>r²</th><th>p-value</th></tr></thead><tbody><tr><td>Bathymetry</td><td>Mean depth [log; m]</td><td>0.28</td><td>0.01</td></tr><tr><td>Photic ratio</td><td>phot<sub>ratio</sub> = z<sub>eul</sub>/z<sub>max</sub> [unitless]</td><td>0.59</td><td>6 × 10<sup>-5</sup></td></tr><tr><td rowspan="3">Organic matter supply</td><td>Drainage ratio [log; unitless]</td><td>0.03</td><td>0.49</td></tr><tr><td>NDVI [log; unitless]</td><td>0.12</td><td>0.12</td></tr><tr><td>Macrophyte coverage [%]</td><td>0.22</td><td>0.07</td></tr><tr><td>Sediment heat release</td><td>Warming rate under the ice, K<sub>heat</sub> [log; d<sup>-1</sup>]</td><td>0.33</td><td>8 × 10<sup>-3</sup></td></tr></tbody></table>

indicated a strong relationship tying the DO decay rate K<sub>exp</sub> to the deep warming rate K<sub>heat</sub>, thus pointing, for cryostratified lakes, to a key role of the sediment heat release on the DO dynamics under ice. A testbed model explored whether differences in the near-bed stratification could be the mechanism by which K<sub>exp</sub> depended on varying deep warming rates under the ice but also potentially on the lakes cryoregimes.

## Modeling the DO dynamics under the ice

Near-bottom DO decreased after ice-on, with a high depletion rate over the first days and weakening as the ice-on period progressed, consistently with the observations by Terzhevik et al. (2009) and Smits et al. (2021). The high-resolution monitoring of near-bottom DO also revealed episodes of complex patterns and deep oxygen resupply that caused DO evolution to diverge from an exponential decay. On two occasions (in Bresses inférieur and Merlet), the disruptions in DO patterns forced us to limit the fit of the exponential model to the period prior to DO disruption. Investigating whether disruptions in the exponential DO decay arise from technical aspects (vertical shift of the mooring position) or physical processes (e.g., under-ice convection, snowmelt, or ice or snow subsidence) would require additional data such as light, pressure and weather data (Smits et al. 2021), that were not available in this study. However, DO decay resumed quickly after temporary episodes of oxygen resupply, and the exponential DO decay model, even when computed only from a part of the dataset, fitted well with the observation throughout the iced season for all but one lake: for Bresses Inferieur, the temporal trajectory deviated clearly from an exponential curve after about 50 d. The time to reach hypoxia lagged by 100 d compared to what could be predicted from K<sub>exp</sub> computed over the first 50 d. The disruption in DO dynamics in Lake Bresses Inférieur coincided with a complex thermal dynamic of the water column from day 50 due to the intrusion of colder water at the lake bottom (SI 1), pointing to a resupply of oxygen likely due to ground-water inflows rather than episodic surface events.

The predominant DO decay dynamics in all lakes and years throughout winter confirms that the Winter stage I sensu Kirillin et al. (2012), i.e., snow-covered ice, is the most

common model for the iced season of high-altitude mountain lakes, echoing conclusions from Smits et al. (2021). In such conditions, K<sub>exp</sub>, the exponential DO decay rate, appears as a handy universal metric to quantify, compare and model the near-bottom DO dynamics throughout the ice-cover season. The duration of hypoxia under ice was much more strongly dependent on DO decay rate K<sub>exp</sub> than on ice duration or initial DO concentrations, and the time to reach hypoxia was accurately predicted from K<sub>exp</sub>. Addressing the drivers of K<sub>exp</sub> variability among and within lakes and the processes governing DO decay is key to predicting under-ice near-bottom hypoxia for alpine lakes, for current conditions and under future climate change.

## Physical controls of DO decay under the ice

K<sub>exp</sub> varied by more than one order of magnitude among lakes and by a factor 1.8 among years for the lake represented by more than 2 yrs in the dataset (Lake Anterne). In our study, linear rates of change computed in the first weeks after the ice onset were, on average, equivalent to 0.11 mg L<sup>-1</sup> d<sup>-1</sup>, except at Lake Mont Coua, where it was 0.30 mg L<sup>-1</sup> d<sup>-1</sup>. These values are well within the reported range of 0.10–0.40 mg L<sup>-1</sup> d<sup>-1</sup> found at the bottom of other high-altitude lakes (Granados et al. 2020; Smits et al. 2021) or low-land oligotrophic lakes (Ghane and Boegman 2023). The high spatial variability in K<sub>exp</sub> was less surprising than its temporal variation considering the already well-observed influence of lake morphology on the DO dynamics under ice (Mathias and Barica 1980; Leppi et al. 2016). The DO decay rates are typically faster in shallower lakes under ice (Mathias and Barica 1980; Leppi et al. 2016) as a result of a greater contribution of sediment processes to the water column conditions, i.e., greater sediment surface to water volume ratio (Bengtsson and Ali-Maher 2020).

Sediment oxygen uptake, which is ultimately the process driving DO consumption under ice, is expected to increase along with the sediment OM content and therefore with lake productivity (Steinsberger et al. 2017). However, surprisingly, the values of K<sub>exp</sub> in the poorly productive lakes of this study were of the same order of magnitude as those observed in the mesotrophic, organic-rich Lake Vendyurskoe in Russia (Terzhevik et al. 2009). Furthermore, K<sub>exp</sub> did not vary along with estimators for OM inputs to the sediment. Although proper data on sediment OM would be necessary, our study points to SOU rates having a second role in explaining the variability in K<sub>exp</sub> within and among the study lakes. The variability of the near-bed DO dynamics that we observed among and within lakes primarily arises from physical processes rather than changes in biogeochemical rates. Our results point to a crucial role of sediment heat release, and more broadly to the key role played by the near-bed stratification in modulating the rate of DO decay under the ice of low-productivity mountain lakes.
<figure><img src="image_9.png" /><figcaption>**Fig. 7.** (a) Effect of the cryoregime and sediment heat flux at ice-on (i.e., release) on *K*<sub>exp</sub> at 1 m above the sediment for a 20 m deep lake. (b–c) Changes in the deep water temperature (14 m depth to the bottom) from the ice onset to the first month of ice cover for a cryomictic (b) and cryostratified (c) lake for a null sediment heat flux (purple) or a sediment heat flux at ice onset of 4 W m<sup>-2</sup> (yellow). (d and e) Changes in the deep water DO for a cryomictic (d) and cryostratified (e) from the ice onset to the first month of ice cover for a null sediment heat flux (purple) or a sediment heat flux at ice onset of 4 W m<sup>-2</sup> (yellow).</figcaption></figure>

Linking near-bed stratification to DO decay through a
testbed model

In cryostratified lakes, DO decay was significantly faster in
lakes with higher deep warming rates. For the cryostratified

Lake Anterne, for which we had 4 yrs of monitoring, *K*<sub>exp</sub> var-
ied significantly among years, where *K*<sub>exp</sub> increased 1.8-fold for
a tripling of *K*<sub>heat</sub>. Our results align with the observations of
Mortimer and Mackereth (1958) and simulations by Golosov
et al. (2007) and Terzhevik et al. (2009), theorizing a mecha-
nistic link between the deep water warming from sediment
heat release and DO dynamics under ice. Using a reaction–
diffusion model, we provided a proof of concept that near-bed
stratification (i.e., near-bottom density gradient) generated by
sediment heat release can produce different DO decay rates by
limiting DO diffusion from upper water layers.

Our proof-of-concept does not exclude other mechanisms
that have been previously proposed, but not yet formally
proven, to underly the relationships between DO decay rate
and sediment heat. For instance, we did not account for heat
or salts advection from downslope flowing currents (Golosov
et al. 2007; Terzhevik et al. 2009; MacIntyre et al. 2018), all of
which should nevertheless reinforce the near-bottom density
gradient. Downslope currents could also complexify the pic-
ture if they result in significant cross-shore advection, which
has been observed but not yet quantified (Stefanovic and Ste-
fan 2002; Ramón et al. 2021). Nor did we include the
suggested temperature-dependency of DO gradient within the
sediment (Golosov et al. 2007) due to the lack of observations
to build upon. Although greatly simplified, the model pro-
vides evidence that the strength of stratification near the lake
bottom, resulting from varying heat sediment release among
lakes or years, can significantly contribute to the near-bottom
oxygen dynamics by modulating DO vertical diffusion.

The K<sub>exp</sub> observed for the cryomictic lake stood out from
the whole dataset, reaching a value that was twice as large as
those observed for cryostratified lakes of similar depth. Our
simulation results also suggest that near-bed stratification
could explain the greater DO decay rates observed in
cryomictic lakes compared to all other cryostratified lakes.
Cryomictic and cryostratified lakes differ in their depth-
averaged ice-on temperatures and strength of the inverse strat-
ification under the ice (Yang et al. 2021). For cryostratified
lakes, low wind exposure in fall restricts the heat loss to the
water surface, until the surface reaches the congelation point
(Yang et al. 2021). The temperature and thereby density gradi-
ent of the inverse stratification is located close to the surface.
In contrast, cryomictic lakes are exposed to wind mixing so
that a much thicker water layer has to lose heat before the sur-
face reaches the freezing point (Yang et al. 2021) and the
temperature gradient is much closer to the bottom. A strong
near-bed temperature gradient occurs in cryomictic lakes
regardless of sediment heat release, with strong consequences
on the DO diffusive fluxes and thereby K<sub>exp</sub>. The proposed phys-
ical model suggests a critical role of lake cryoregimes on under-
ice DO dynamics, calling for observational field confirmation.

Physical controls of DO decay under ice and the light
climate of open waters

Our results indicate that the summer light regime is a criti-
cal predictor of oxygen dynamics under the ice of
cryostratified lakes. Specifically, the photic ratio (measured in
the summer season) was the best predictor of oxygen decay

under ice. As sediment heat release is important for DO decay
rates, the factors controlling sediment heat accumulation dur-
ing open waters, such as light climate, indirectly control DO
dynamics under ice (Golosov et al. 2007). Sediment heat accu-
mulation depends on water transparency (Fang and Ste-
fan 1996) and lake hypsometry, constraining the area of
sediment directly exposed to solar irradiance. The photic ratio,
measured at the end of summer, varies significantly among
lakes and years, with photic depths as shallow as 5% of the
maximum lake depth for lakes receiving glacial flour (Arpont,
Muzelle for 2018–2019) and as deep as 220% of the maximum
depth for the clearest ones (Noir du Carro, Pisses). Light pen-
etration down to the bottom during summer generates weak
density stratification with a greater heat accumulation in the
sediment. In contrast, sediment heat accumulation will be
restricted to the shallower sediments in more turbid lakes
when more stable stratified waters prevent downward heat
penetration (Fang and Stefan 1996). Lakes of greater photic
depth might have accumulated considerable amounts of heat
in summer, and indeed we found a positive correlation
between the lake photic ratio and K<sub>heat</sub> (r = 0.44, p = 0.04).
Greater sediment heat storage in summer generates a stronger
near-bed density gradient when heat is released under the ice
(Golosov et al. 2007), causing, ultimately, a faster DO decay at
the lake bottom during the ice-covered period.

DO under the ice of alpine lakes and climate change

In our study lakes, both the time to reach hypoxia and the
duration of hypoxia under ice could be predicted well from
the DO decay rate itself varying primarily, for cryostratified
lakes, with the photic ratio. Ultimately, the combined lake
depth and water transparency could also very well explain the

<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Hypoxia duration (in days) as a function of the photic ratio (open symbol = cryomictic, closed symbols = cryostratified). The linear fit is performed for cryostratified lakes only. The equation of the linear model is Hypoxia_duration = 69.7 Photic_ratio + 20.5. The size of data points is proportional to the weight given within the linear model.</figcaption></figure>
variability in the duration of hypoxia among and within lakes (Fig. 8), providing a framework through which the effect of climate change on the under-ice DO dynamics of alpine cryostratified lakes can be theorized.

While reduced in duration, the ice cover of lakes in the European Alps at altitude above 1500 m asl is predicted to persist until the end of the century even in the worst-case climate scenarios (Råman Vinnå et al. 2021). We expect changes in the duration of hypoxia under the ice to be more drastic in the shallowest lakes, although the magnitude might be case-specific. While the water transparency of mountain lakes is expected to undergo significant changes under climate warming (Parker et al. 2008), those changes are not expected to be monotonic, but rather dependent on the specific landscape features (Olson et al. 2017). Water transparency in mountain lakes responds mainly to catchment inputs in turbidity (for glacially-fed lakes) and colored dissolved organic matter (CDOM, for non-glacially-fed lakes). With glacier retreat, lakes get progressively disconnected from glacial inputs, leading to a turbidity of glacial-fed lakes that decreases and gets more variable in time (Olson et al. 2017; Olson et al. 2020). As disconnection gets complete, lakes suddenly shift from turbid to clearwater with substantial consequences on the heat content and vertical distribution in summer (Perga et al. 2018). Assuming a stable cryostratified regime, and a Secchi depth of 0.5 m for turbid glacial-fed lakes (as of Arpont and Muzelle 2018–2019) and 5 m for a non-glacially fed lake, a full glacial disconnection shall lead to hypoxia lasting 53 d longer in a 10 m max depth lake, and 21 d longer in a 30 m depth lake (based on Fig. 8). In both cases, this predicted lengthening of hypoxia is of similar scale and may outcompete the predicted reduction of ice cover for high altitude lakes of the Alps by the end of the century for the RCP2.6 scenario (20 d, e.g., Råman Vinnå et al. 2021). Hypoxia duration under the ice might remain steady or potentially increase under the ice of lakes getting disconnected from glacial waters in the future.

For alpine lakes that are currently non-glacially-fed, or for those that will lose glacial connection, CDOM controls diffuse light attenuation (Olson et al. 2017). Catchment greening (Carlson et al. 2017) is expected to promote CDOM exports to lakes (Olson et al. 2017). However, CDOM absorbs mainly in the UV domain, so that PAR attenuation does not scale to CDOM concentrations in the range of values observed for mountain lakes (Olson et al. 2017). On top, neither the lake Secchi depth nor the photic ratio were related to the catchment NDVI in our dataset, highlighting that the relationship between light climate and catchment productivity is far from direct in alpine lakes. Increased CDOM loads to mountain lakes might therefore not significantly affect their photic ratio. Sediment heat accumulation might nevertheless increase under a warmer climate, given constant water clarity. Consequences of climate change on the under-ice bottom DO dynamics in alpine lakes will then depend on which of the ice

cover duration or the structure of the inverse stratification takes over.

## Data Availability Statement

The sensors data that support the findings of this study are openly available in IS OLA (Observatory on LAkes) at https://si-ola.inrae.fr/si_lacs/login.jsf. Processed data are available at https://doi.org/10.5281/zenodo.7674243.

## References

Bengtsson, L., and O. Ali-Maher. 2020. The dependence of the consumption of dissolved oxygen on lake morphology in ice covered lakes. Hydrol. Res. **51**: 381–391. doi:10.2166/ nh.2020.150

Beutel, M. W. 2003. Hypolimnetic anoxia and sediment oxygen demand in California drinking water reservoirs. Lake Reserv. Manag. **19**: 208–221. doi:10.1080/07438140309354086

Birck, C., and others. 2013. Sentinel lakes: A network for the study and management of mountain lakes in the French Alps and in Corsica. Eco Mont J. Prot. Mt. Areas Res. **5**: 63–69. doi:10.1553/eco.mont-5-1s63

Bouffard, D., J. D. Ackerman, and L. Boegman. 2013. Factors affecting the development and dynamics of hypoxia in a large shallow stratified lake: Hourly to seasonal patterns. Water Resour. Res. **49**: 2380–2394. doi:10.1002/wrcr.20241

Brown, L. C., and C. R. Duguay. 2010. The response and role of ice cover in lake-climate interactions. Prog. Phys. Geogr. Earth Environ. **34**: 671–704. doi:10.1177/0309133310375653

Caldwell, T. J., and others. 2021. Drivers and projections of ice phenology in mountain lakes in the western United States. Limnol. Oceanogr. **66**: 995–1008. doi:10.1002/lno.11656

Carlson, B. Z., M. C. Corona, C. Dentant, R. Bonet, W. Thuiler, and P. Choler. 2017. Observed long-term greening of alpine vegetation—a case study in the French Alps. Environ. Res. Lett. **12**: 114006. doi:10.1088/1748-9326/aa84bd

Cavaliere, E., and others. 2021. The Lake ice continuum concept: Influence of winter conditions on energy and ecosystem dynamics. Journal of geophysical research. Biogeosciences **126**: e2020JG006165. doi:10.1029/2020JG006165

Couture, R.-M., H. A. de Wit, K. Tominaga, P. Kiuru, and I. Markelov. 2015. Oxygen dynamics in a boreal lake responds to long-term changes in climate, ice phenology, and DOC inputs. J. Geophys. Res. Biogeo. **120**: 2441–2456. doi:10.1002/2015JG003065

Cremona, F., A. Laas, P. C. Hanson, M. Sepp, P. Nõges, and T. Nõges. 2019. Drainage ratio as a strong predictor of allochthonous carbon budget in Hemiboreal Lakes. Ecosystems **22**: 805–817. doi:10.1007/s10021-018-0304-7

Davis, M. N., T. E. McMahon, K. A. Cutting, and M. E. Jaeger. 2020. Environmental and climatic factors affecting winter hypoxia in a freshwater lake: Evidence for a hypoxia refuge and for re-oxygenation prior to spring ice loss.
Hydrobiologia **847**: 3983-3997. doi:10.1007/s10750-020-04382-z

Devlin, S. P., M. J. Vander Zanden, and Y. Vadeboncoeur. 2016. Littoral-benthic primary production estimates: Sensitivity to simplifications with respect to periphyton productivity and basin morphometry. Limnol. Oceanogr.: Methods **14**: 138-149. doi:10.1002/lom3.10080

Fang, X., and H. G. Stefan. 1996. Dynamics of heat exchange between sediment and water in a lake. Water Resour. Res. **32**: 1719-1727. doi:10.1029/96WR00274

Flaim, G., D. Andreis, S. Piccolroaz, and U. Obertegger. 2020. Ice cover and extreme events determine dissolved oxygen in a Placid Mountain Lake. Water Resour. Res. **56**: e2020WR027321. doi:10.1029/2020WR027321

Fountain, A. G., J. L. Campbell, E. A. G. Schuur, S. E. Stammerjohn, M. W. Williams, and H. W. Ducklow. 2012. The disappearing cryosphere: Impacts and ecosystem responses to rapid cryosphere loss. BioScience **62**: 405-415. doi:10.1525/bio.2012.62.4.11

Gammons, C. H., and others. 2014. Stable isotopes track bio-
geochemical processes under seasonal ice cover in a shal-
low, productive lake. Biogeochemistry **120**: 359-379. doi:
10.1007/s10533-014-0005-z

Ghane, A., and L. Boegman. 2023. The dissolved oxygen bud-
get of a small Canadian shield lake during winter. Limnol.
Oceanogr. **68**: 265-283. doi:10.1002/lno.12265

Gobiet, A., S. Kotlarski, M. Beniston, G. Heinrich, J. Rajczak, and M. Stoffel. 2014. 21st century climate change in the European Alps—a review. Sci. Total Environ. **493**: 1138–1151. doi:10.1016/j.scitotenv.2013.07.050

Golosov, S., and G. Kirillin. 2010. A parameterized model of heat storage by lake sediments. Environ. Model. Softw. **25**: 793-801. doi:10.1016/j.envsoft.2010.01.002

Golosov, S., O. A. Maher, E. Schipunova, A. Terzhevik, G. Zdorovennova, and G. Kirillin. 2007. Physical background of the development of oxygen depletion in ice-covered lakes. Oecologia **151**: 331-340. doi:10.1007/s00442-006-0543-8

Granados, I., M. Toro, S. Giralt, A. Camacho, and C. Montes. 2020. Water column changes under ice during different winters in a mid-latitude Mediterranean high mountain lake. Aquat. Sci. **82**: 30. doi:10.1007/s00027-020-0699-z

Hampton, S. E., and others. 2017. Ecology under lake ice. Ecol. Lett. **20**: 98-111. doi:10.1111/ele.12699

Heinz, G., J. Ilmberger, and M. Schimmele. 1990. Vertical mixing in Überlinger see, western part of Lake Constance. Aquat. Sci. **52**: 256-268. doi:10.1007/BF00877283

Hondzo, M., and G. Stefan. 1993. Lake water temperature sim-
ulation model. J. Hydraul. Eng. **119**: 1251-1273. doi:10.
1061/(ASCE)0733-9429(1993)119:11(1251)

Jansen, J., and others. 2021. Winter limnology: How do
hydrodynamics and biogeochemistry shape ecosystems
under ice? Journal of geophysical research. Biogeosciences
**126**: e2020JG006237. doi:10.1029/2020JG006237

Jassby, A., and T. Powell. 1975. Vertical patterns of eddy diffu-
sion during stratification in castle Lake, California. Limnol.
Oceanogr. **20**: 530-543. doi:10.4319/lo.1975.20.4.0530

Joint Research Centre, Institute for Environment Sustainability. 2012. *In* S. Poikane [ed.], Water framework directive inter-calibration technical report. Part 2: Lakes. Publications Office. 10.2788/23415

Kincaid, D. W., E. C. Adair, D. Joung, J. D. Stockwell, and A. W. Schroth. 2022. Ice cover and thaw events influence nitrogen partitioning and concentration in two shallow eutrophic lakes. Biogeochemistry **157**: 15-29. doi:10.1007/s10533-021-00872-x

Kirillin, G., and others. 2012. Physics of seasonally ice-covered lakes: A review. Aquat. Sci. **74**: 659-682. doi:10.1007/s00027-012-0279-y

Larsen, S., T. Andersen, and D. O. Hessen. 2011. Predicting organic carbon in lakes from climate drivers and catchment properties. Global Biogeochem. Cycles **25**: GB3007. doi:10.1029/2010GB003908

Leppi, J. C., C. D. Arp, and M. S. Whitman. 2016. Predicting late winter dissolved oxygen levels in Arctic Lakes using morphology and landscape metrics. Environ. Manag. **57**: 463-473. doi:10.1007/s00267-015-0622-x

MacIntyre, S., A. Cortés, and S. Sadro. 2018. Sediment respi-
tion drives circulation and production of CO₂ in ice-
covered Alaskan arctic lakes. Limnol. Oceanogr. Lett. **3**:
302-310. doi:10.1002/lo12.10083

Magnuson, J. J., A. L. Beckel, K. Mills, and S. B. Brandt. 1985. Surviving winter hypoxia: Behavioral adaptations of fishes in a northern Wisconsin winterkill lake. Environ. Biol. Fishes **14**: 241-250. doi:10.1007/BF00002627

Magnuson, J. J., and others. 2000. Historical trends in lake and river ice cover in the northern hemisphere. Science **289**: 1743-1746. doi:10.1126/science.289.5485.1743

Malm, J., and others. 1998. Field study on currents in a shal-
low, ice-covered lake. Limnol. Oceanogr. **43**: 1669-1679.
doi:10.4319/lo.1998.43.7.1669

Mathias, J. A., and J. Barica. 1980. Factors controlling oxygen depletion in ice-covered lakes. Can. J. Fish. Aquat. Sci. **37**: 185-194. doi:10.1139/f80-024

Mortimer, C. H., and F. J. H. Mackereth. 1958. Convection and its consequences in ice-covered lakes, p. 923-932. In SIL Proceedings, 1922-2010, Taylor and Francis online; **13**: 923-932. doi:10.1080/03680770.1956.11895490

Müller, B., L. D. Bryant, A. Matzinger, and A. Wüest. 2012. Hypolimnetic oxygen depletion in eutrophic lakes. Environ. Sci. Technol. **46**: 9964-9971. doi:10.1021/es301422r

Obertegger, U., B. Obrador, and G. Flaim. 2017. Dissolved oxygen dynamics under ice: Three winters of high-frequency data from Lake Tovel, Italy. Water Resour. Res. **53**: 7234-7246. doi:10.1002/2017WR020599

Olson, M. H., J. M. Fischer, C. E. Williamson, E. P. Overholt, and N. Theodore. 2017. Landscape-scale regulators of water transparency in mountain lakes: Implications of projected
glacial loss. Can. J. Fish. Aquat. Sci. **75**: 1169-1176. doi: **10.1139/cjfas-2017-0215**

Olson, M. H., J. M. Fischer, M. Hayashi, and C. E. Williamson. 2020. Meteorological drivers of interannual variation in transparency of mountain lakes. Arct. Antarct. Alp. Res. **52**: 424-434. doi:10.1080/15230430.2020.1800972

Parker, B. R., R. D. Vinebrooke, and D. W. Schindler. 2008. Recent climate extremes alter alpine lake ecosystems. Proc. Natl. Acad. Sci. USA **105**: 12927-12931. doi:10.1073/pnas.0806481105

Perga, M. E., R. Bruel, L. Rodriguez, Y. Guenand, and D. Bouffard. 2018. Storm impacts on alpine lakes: Antecedent weather conditions matter more than the event intensity. Glob. Chang. Biol. **24**: 5004-5016. doi:10.1111/gcb.14384

Preston, D. L., and others. 2016. Climate regulates alpine lake ice cover phenology and aquatic ecosystem structure. Geophys. Res. Lett. **43**: 5353-5360. doi:10.1002/2016GL069036

Råman Vinnå, L., I. Medhaug, M. Schmid, and D. Bouffard. 2021. The vulnerability of lakes to climate change along an altitudinal gradient. Commun. Earth Environ. **2**: 35. doi: **10.1038/s43247-021-00106-w**

Ramón, C. L., H. N. Ulloa, T. Doda, K. B. Winters, and D. Bouffard. 2021. Bathymetry and latitude modify lake warming under ice. Hydrol. Earth Syst. Sci. **25**: 1813-1825. doi:10.5194/hess-25-1813-2021

Salonen, K., M. Leppäranta, M. Viljanen, and R. D. Gulati. 2009. Perspectives in winter limnology: Closing the annual cycle of freezing lakes. Aquat. Ecol. **43**: 609-616. doi:10.1007/s10452-009-9278-z

Sharma, S., and others. 2021. Loss of ice cover, shifting phe-nology, and more extreme events in northern hemisphere lakes. J. Geophys. Res.: Biogeosci. **126**: e2021JG006348. doi:10.1029/2021JG006348

Smits, A. P., N. W. Gomez, J. Dozier, and S. Sadro. 2021. Win-ter climate and lake morphology control ice phenology and under-ice temperature and oxygen regimes in mountain lakes. J. Geophys. Res.: Biogeosci. **126**: e2021JG006277. doi:10.1029/2021JG006277

Song, S., and others. 2019. Under-ice metabolism in a shallow lake in a cold and arid climate. Freshw. Biol. **64**: 1710-1720. doi:10.1111/fwb.13363

Stefanovic, D. L., and H. G. Stefan. 2002. Two-dimensional temperature and dissolved oxygen dynamics in the littoral region of an ice-covered lake. Cold Reg. Sci. Technol. **34**: 159-178. doi:10.1016/S0165-232X(02)00003-4

Steinsberger, T., M. Schmid, A. Wüest, R. Schwefel, B. Wehrli, and B. Müller. 2017. Organic carbon mass accumulation rate regulates the flux of reduced substances from the sedi-ments of deep lakes. Biogeosciences **14**: 3275-3285. doi:10.5194/bg-14-3275-2017

Stepanenko, V., I. Mammarella, A. Ojala, H. Miettinen, V. Lykosov, and T. Vesala. 2016. LAKE 2.0: A model for tem-perature, methane, carbon dioxide and oxygen dynamics in lakes. Geosci. Model Dev. **9**: 1977-2006. doi:10.5194/gmd-9-1977-2016

Terry, J. A., A. Sadeghian, and K.-E. Lindenschmidt. 2017. Modelling dissolved oxygen/sediment oxygen demand under ice in a shallow eutrophic prairie reservoir. Water **9**: 131. doi:10.3390/w9020131

Terzhevik, A., and others. 2009. Some features of the thermal and dissolved oxygen structure in boreal, shallow ice-covered Lake Vendyurskoe, Russia. Aquat. Ecol. **43**: 617-627. doi:10.1007/s10452-009-9288-x

Ulloa, H. N., K. B. Winters, A. Wüest, and D. Bouffard. 2019. Differential heating drives downslope flows that accelerate mixed-layer warming in ice-covered waters. Geophys. Res. Lett. **46**: 13872-13882. doi:10.1029/2019GL085258

Wang, X., and others. 2021. High-resolution mapping of ice cover changes in over 33,000 lakes across the north temper-ate zone. Geophys. Res. Lett. **48**: e2021GL095614. doi:10.1029/2021GL095614

Woolway, R. I., and C. J. Merchant. 2019. Worldwide alter-ation of lake mixing regimes in response to climate change. Nat. Geosci. **12**: 271-276. doi:10.1038/s41561-019-0322-x

Yang, B., and others. 2021. A new thermal categorization of ice-covered lakes. Geophys. Res. Lett. **48**: e2020GL091374. doi:10.1029/2020GL091374

Zdorovennova, G., and others. 2021. Dissolved oxygen in a shallow ice-covered lake in winter: Effect of changes in light, thermal and ice regimes. Water **13**: 2435. doi:10.3390/w13172435

## Acknowledgments

We thank all three reviewers for their constructive comments on our manuscript. None of this would have been possible without the involve-ment of the managers and rangers of the "Parcs nationaux et espaces protégés de France" and the extraordinary energy of the "reseau lacs sen-tilinelles." We would like to thank the scientific managers of ASTERS (C. Birk), the Parc National de la Vanoise (V. Augé), the Parc National des Ecrins (R. Bonnet and C. Sagot), the Parc du Mercantour (MF Leccia) and EDF (V. Chanudet), as well as all their involved staff. Open access funding provided by Universite de Lausanne.

## Conflict of Interest

None declared.

Submitted 31 May 2022
Revised 16 December 2022

Accepted 08 March 2023

Associate editor: Leon Boegman