Emerging unprecedented lake ice loss in climate change projections

The study addresses how greenhouse warming alters lake ice phenology and thickness globally, and when anthropogenic signals will emerge from natural variability. Recent decades have seen widespread loss of lake ice, with delayed freeze and earlier breakup. Prior work linked observed changes to greenhouse forcing and projected continued loss but often used offline lake models or statistical approaches that do not fully capture thickness changes, diurnal thermodynamics, or two-way coupling between lakes, ice, and the atmosphere. This study employs an interactive lake simulator within a fully coupled Earth system model and a large ensemble strategy to isolate forced responses from internal variability, aiming to quantify global and regional trends in lake ice phenology and thickness and to determine the timing and warming levels at which no-analog ice conditions emerge for aquatic ecosystems.

Previous research documented widespread reductions in lake ice cover across the Northern Hemisphere and attributed changes to anthropogenic warming using offline lake models driven by climate model outputs. Statistical models based on empirical relationships between ice phenology and meteorology project intermittent or permanent winter ice loss by century’s end but typically lack projections of ice thickness and detailed duration dynamics. Offline, process-based lake-ice models resolve thermodynamics but miss coupled lake–ice–atmosphere feedbacks and diurnal cycles, which are important especially for large lakes that influence local climate. Observations indicate unprecedented recent rates of ice loss and significant ecological impacts, including altered habitats and risks of local extirpation for cold-adapted species under future no-analog ice regimes. The value of initial-condition large ensembles has been highlighted for disentangling forced trends from natural variability, yet previous multimodel ensemble approaches can suppress variability and mix sensitivities across models, complicating emergence assessments.

The authors analyze output from the CESM2 Large Ensemble (CESM2-LE), consisting of 100 initial-condition perturbation members spanning 1850–2100 under CMIP6 historical forcing through 2014 and SSP3-7.0 thereafter. A subset of 90 members with daily lake ice thickness is used to compute the forced signal (90-member ensemble mean) and natural variability (ensemble spread). Lakes in Antarctica and Greenland are excluded. Ice phenology is defined for lakes with more than 5 days of ice cover per year; freeze date is the first day with ice (thickness > 0 m), breakup the last day before ice-free conditions, and duration the total days with ice. Model framework: The lake component is CLM5’s Lake, Ice, Snow, and Sediment Simulator (LISSS), a one-dimensional thermodynamic model with vertically resolved snow (up to 5 layers), lake water (10 layers including liquid and ice fractions), soil (10 layers), and bedrock (5 layers). It employs bulk transfer for surface fluxes and includes wind stirring, convection, molecular diffusion, and unresolved mixing. Ice albedo decreases during melt to represent puddling/disintegration; snow albedo follows non-vegetated land formulations. The default timestep (1800 s) resolves diurnal cycles. Each 0.9°×1.25° grid cell contains one representative lake with fixed volume layers determined by nominal depth; lake depth and area fraction are aggregated from high-resolution global datasets. Validation: Ice phenology is validated against 331 observational records (≥20 years since 1850; analysis emphasized records ≥20 years since 1950) from the Global Lake and River Ice Phenology database. Spatial variance and climatological means are well reproduced: correlations r = 0.94 (ice duration), 0.78 (freeze date), and 0.94 (breakup date); MAE ≈ 17, 12, and 16 days, respectively. Ice thickness is validated using 15 Canadian records (1947–2002), with spatial correlation r = 0.93 for annual maximum thickness, though the model underestimates variance toward thicker ice (regression slope ≈ 2.1). A known warm bias (~1.5 °C in 30°N–70°N during 1981–2020 versus ERA5) contributes to shorter simulated ice duration. Time of emergence and no-analog conditions: For each lake, the natural variability baseline is defined from 1850–1950 as the mean µ and ±2σ range across time and the ensemble. Emergence of forced change (t0) occurs when the ensemble mean minus 2σ falls below the lower 2σ limit of the natural range for at least five consecutive years. No-analog conditions (t1) are defined when the ensemble mean plus 2σ falls below the lower 2σ limit, implying only a 2.5% chance of returning to the natural range under internal variability. Adaptation time is t1 − t0. Statistical analyses (correlations, trends) use ensemble-mean fields; linear trends and changes are computed for 2000–2100 or specified intervals. Additional attribution of drivers employs subsets of 60 members with necessary outputs to diagnose seasonal changes in sea ice, surface air temperature, lake surface albedo, and absorbed shortwave radiation.

- Global trends (2020–2100): Lakes are projected to freeze 20 ± 8 days later and break up 20 ± 7 days earlier, shortening ice duration by 38 ± 11 days on average. Global mean annual maximum lake ice thickness declines by 0.23 ± 0.07 m over the next 80 years, with the forced signal exiting the envelope of natural variability by 2000–2010. - Regional patterns: Strongest shortening of ice duration (> 0.45 days/year loss) and thinning (≈ −0.04 m/decade) occur in the Canadian Arctic, northern Siberian and Barents coasts, and on the Tibetan Plateau. Over Norway and Sweden, reductions are smaller (≈ −0.28 days/year). In the Southern Hemisphere (outside Antarctica), thin, short-lived ice near the southern Andes becomes mostly intermittent or disappears by century’s end. - Intermittent versus perennial ice: An additional 11.7 ± 1.6% of global lake surface area transitions from climatological winter ice coverage to intermittent coverage by late 21st century. - Temperature–ice sensitivity: Simulated anomalies of ice duration correlate strongly with annual mean surface air temperature (r = −0.8, P < 0.001). The model’s global-scale sensitivity is ≈ 9.9 days decrease in ice duration per 1 °C local warming, close to observational sensitivity ≈ 8.9 days/°C (observational correlation r ≈ −0.47). - Drivers of rapid regional loss: Canadian Arctic lakes experience accelerated loss linked to Arctic and Hudson Bay sea-ice retreat and cold-season polar amplification (heat capacitor feedback). Tibetan Plateau lakes show rapid spring–autumn ice loss amplified by strong insolation and a pronounced bi-seasonal ice–albedo feedback, enhancing heat absorption and delaying freeze. - Emergence and no-analog conditions: First anthropogenic deviations from natural ice duration ranges occur in 1980–1990. No-analog lake ice conditions begin to emerge when global mean temperature reaches about +1.9 °C relative to 1850–1950. Many regions, notably the Canadian Arctic and Tibetan Plateau, are projected to encounter no-analog ice coverage within 4–5 decades under SSP3-7.0. - Ecological implications: The projected 38 ± 11 day shortening in ice duration aligns with expected extensions of thermal stratification (≈ +33 ± 12 days in related studies) and is consistent with observed shifts in phytoplankton assemblages and other ecological indicators tied to diminishing ice cover.

By leveraging a large ensemble of a single, fully coupled Earth system model with an interactive lake simulator, the study robustly separates anthropogenic signals from internal variability in lake ice dynamics. It identifies when and where forced trends emerge and cross into no-analog territory for aquatic ecosystems. The analysis highlights two hotspots—Canadian Arctic and Tibetan Plateau—where physical mechanisms accelerate ice loss: cold-season polar amplification mediated by sea-ice retreat in the Arctic, and strong seasonal ice–albedo feedback under high insolation on the Tibetan Plateau. These physical drivers explain the accelerated regional trends and clarify how warming translates into earlier breakup, later freeze, and thinning. The findings have broad ecological relevance, as ice phenology orchestrates seasonal light, temperature, mixing, nutrient cycling, and oxygen dynamics, influencing food webs, bloom timing, and species distributions. The timing of emergence corresponds with observed biotic regime shifts, lending confidence to the projections and underscoring the urgency of adaptation and mitigation strategies to safeguard cold-adapted species and ecosystem services dependent on lake ice.

The study provides a physically consistent, mechanistic projection of global and regional lake ice loss using a 90-member large ensemble with an interactive lake model. It quantifies substantial late freeze, early breakup, shortened duration (−38 ± 11 days), and thinning (−0.23 ± 0.07 m) by 2100, identifies hotspots and their drivers, and constrains the timing and warming levels for no-analog lake ice conditions, which begin near +1.9 °C global warming. These results emphasize significant ecological risks under sustained high emissions and support mitigation to limit global mean temperatures below 1.9 °C to avoid unprecedented lake ice changes. Future research directions include single-forcing ensemble experiments to attribute roles of greenhouse gases versus aerosols, improved representation of lake heterogeneity (multiple lakes per grid, horizontal processes), and refined snow–ice parameterizations to reduce biases.

- Spatial representation: One representative lake per 1° grid cell; intra-cell heterogeneity (elevation, morphology, depth, size) is not resolved, potentially misrepresenting individual lake dynamics, especially for large lakes. - Model dimensionality: The one-dimensional lake model omits horizontal circulation, fetch effects, and ice mechanics that influence freeze/breakup in large lakes, potentially biasing climatological means more than variability. - Snow/precipitation uncertainty: Ice growth is sensitive to snow timing/depth and snowfall uncertainties from the atmospheric model, affecting thickness and albedo evolution. - Warm bias: CESM2 exhibits a ~1.5 °C warm bias over 30°N–70°N (1981–2020), translating to shorter simulated ice duration (~−15 days), and contributing to discrepancies with observations. - Scenario and model dependence: Projections are based on a single model and primarily one high-emissions scenario (SSP3-7.0), which may limit scenario breadth; however, emergence assessments are framed versus warming levels to enhance generality. - Validation constraints: Observational comparisons may be affected by shoreline-centric measurements and lake ice heterogeneity; the model tends to underestimate variance toward thicker ice.