Topographically-controlled contribution of avalanches to glacier mass balance in the 21st century

Domain: physics
Published: 2025-11-26
Editorial outcome: minor_revision

Abstract
Glaciers are often located in steep mountain settings and avalanches from surrounding slopes can strongly influence snow accumulation patterns on their surface. This effect has however never been quantified for more than a few glaciers and the impact on the future evolution of glaciers is unclear. We coupled an avalanche and a glacier model to estimate the contribution of avalanches to the accumulation of all glaciers in the world and how this affects their evolution throughout the 21st century. Globally, 3% of the snow accumulation on glaciers comes from avalanches and 1% is removed by avalanches. This net contribution varies between regions and glaciers, with a maximum of 15% for New Zealand. Accounting for avalanches modifies the altitudinal pattern of glacier mass balance and the projected evolution of individual glaciers. The main effects include (1) a longer persistence of small glaciers, with for example three times more ice retained by glaciers smaller than 1 km2 in Central Europe under a low-emissions scenario, and (2) an increased sensitivity of high-elevation accumulation zones to future warming. We anticipate the relative influence of avalanches to increase in the future and advocate for a better monitoring of this process and representation in glacier models.

Introduction
Glaciers are important components of the world’s ‘water towers’ and are currently undergoing major changes due to atmospheric warming. Future projections show that they will lose considerable mass regionally in the coming decades, thereby increasing the vulnerability of downstream regions to future atmospheric warming. However, these projections remain uncertain for individual glaciers due to a lack of understanding and quantification of key processes influencing glacier mass balance, such as supraglacial debris, high-elevation precipitation and avalanching.

Mountainglaciersusuallygainmassviasolidprecipitationfallingintheiraccumulationareathatisthenadvecteddownstreamwithglacierflow.Thus,themassbalanceistraditionallyexpectedtoincreasewithelevation.Forcatchmentsinsteepmountainrangeshowever,therecanbelargemassinputsfrommountainheadwallsatlocalizedportionsoftheglacier,bothintheaccumulationandablationzones,whichleadstoanon-linearelevationdistributionofmassbalance.Theseinputs,whichvaryinsizeandoriginatefromthegravitationalredistributionofsnoworice,cancontributetothepersistenceofglaciersatlowaltitudesandthereforebufferthedepletionofmountainwaterresources.Inthefewlocationswhereon-glacieravalanchinghasbeenestimated,itisoneofthemainsourcesofmassinputs,withestimatedavalanchecontributionsrangingbetween9%and92%oftheannualaccumulation.However,thispresentandfutureinfluencehasneverbeenquantifiedformorethanafewglaciers,withlargeuncertainties,duetoalackoffieldmeasurementsofmasscontributionfromavalanchingandtailoredremotesensingapproaches.

Globalandregional-scalemodelsusuallyexpresstheglaciermassbalanceasafunctionoftemperatureandprecipitation,withanempiricalcorrectionfactorthatscalesprecipitationinputsandattemptstoaccountforunobservedprocessessuchasavalanching,snowredistributionbywindandorographicprecipitation.Thisprecipitationcorrectionfactor,togetherwithatemperaturebiasparameterandthemeltfactorusedinglobal-scaletemperature-indexmodels,needstobecalibratedonaglacier-by-glacierbasis.Inrecentyears,therehasbeenanemergenceofhigh-resolutionlarge-scalegeodeticestimatescharacterizingthestateoftheworld’sglaciers,pavingthewayforamuchmoreaccuratecalibrationofmassbalanceparametersforeachindividualglacierinregionaltoglobal-scalemodels.However,evenwithbettercalibrationdata,thenon-linearityofthemassbalancepatternsfromavalanches(amongotherprocesses; Fig. 1) is impossible to match by these simple models, and leads to an inaccurate representation of the spatial variability of mass balance and, as a result, glacier flow. In such cases the precipitation correction factor can end up being particularly and unrealistically high to account for the excess mass input. A more explicit representation of these non-linear mass balance contributions from avalanching is likely to help constrain the glacier-wide mass balance parameters and therefore improve the simulations from these glacier models, as recently shown for other processes, notably calving or sub-debris melt.

In this study we couple a global glacier model with a gravitational snow redistribution model to estimate the current contribution of avalanches to glacier accumulation for all glaciers in the world. We use this modelling framework to derive the influence of avalanches on the future evolution of individual glaciers and mountain regions throughout the 21st century and analyse this influence relative to model runs that do not explicitly account for avalanches. Finally, we assess how this avalanche contribution will evolve as glaciers retreat in a warming climate.

Results
Substantial contribution of avalanches to glacier accumulation
We use a simple gravitational snow redistribution model run with W5E5v2.0 mean monthly climate inputs over the glacier model calibration period 01/2000-12/2019 to quantify the contribution of avalanches to the annual glacier accumulation (Methods). The snow redistribution model is run at a monthly timescale and removes the snow from areas where the snowpack exceeds a certain threshold determined by the slope. The maximum allowable snow depth is set to decrease exponentially as the slope increases. This implies that our definition of avalanches encompasses any type of gravitational snow redistribution from steep terrain onto lower angle slopes, which includes all magnitudes of events, from drifting snow along steep headwalls to large slab avalanches. We note that this model does not represent wind-driven snow redistribution, which is a particularly difficult process to constrain given its dependence on an accurate representation of wind fields on a glacier, and despite recent advances in quantifying and modelling this effect at relatively small spatial scales.

The glacier-wide avalanche contribution to accumulation can vary considerably, with values above 50%, indicating that more than half of the annual snow accumulation comes from avalanches, and values lower than −50%, indicating that more than half of the snow accumulated on the glacier surface has been removed by avalanches (Fig. 2; Supplementary Fig. 1; Supplementary Table 1). For the glacierized regions of the world described in the Randolph Glacier Inventory v6.0, the regional avalanche contribution is generally positive, with higher values for regions prone to steeper relief and higher precipitation. New Zealand is the region with the highest net avalanche contribution, with a region-wide contribution of 14.6% or 0.49 m w.eq. yr−1, followed by South-East Asia (11.2%; 0.23 m w.eq. yr−1) and Central Europe (9.1%; 0.22 m w.eq. yr−1). The regions with an avalanche contribution close to 0%, the Arctic Canada North and South, Svalbard, Iceland, the Russian Arctic and the Subantarctic and Antarctic Islands, are characterized by relatively flat glaciers and ice caps with few avalanche-prone slopes. At the global scale, we estimate that avalanches contribute 2.0% of the snow accumulation on glaciers, corresponding to 0.03 m w.eq. yr−1, and 20% of the world’s glaciers have more than 5.9% of their accumulation coming from avalanches. In total, avalanches add on average ~29,000 Mt yr−1 of snow on glaciers and remove ~11,000 Mt yr−1. Alaska and South Asia East have the highest total positive and negative on-glacier avalanche contributions of all glacierized regions on Earth, together accounting for more than half of the global positive and nearly one third of the global negative avalanche contributions (Fig. 2b-c; Supplementary Table 1). The Low Latitudes and Southern Andes regions are the only regions with a negative regional avalanche contribution (−2.3%; −0.04 m w.eq. yr−1 and −1.5%; −0.02 m w.eq. yr−1), and 18% of the snow removed by avalanches from glaciers comes from the Low Latitudes despite the relatively small glacier area of the region. These numbers highlight the importance of this gravitational snow removal process for glaciers in the Andes.

The avalanche contribution is smaller for large glaciers (Supplementary Fig. 1), as avalanches tend to only affect the mass balance at the margins of extensive and relatively flat glacier tongues. Glaciers with steeper average slopes tend to be more influenced by avalanches, whether by snow addition or removal, while glaciers with a positive avalanche contribution tend to have a higher percentage of debris cover (Supplementary Figs. 2–4). This confirms past observations conducted at the scale of High Mountain Asia, which likely indicate that debris-covered and avalanche-fed glaciers exist in a similar topographic setting, with steep surrounding slopes. The glacier-wide net avalanche contribution is also negatively correlated with the glacier mean elevation (Supplementary Figs. 2–4). Therefore, avalanches appear to be contributing to maintaining glaciers at lower elevations for a given climate (Supplementary Fig 2-4).

Increased spatial variability of glacier mass balance with avalanches
Glacier models implicitly account for avalanches because they are calibrated against geodetic mass balance observations, which integrate the effects of all the different processes affecting the mass balance. However, the spatial variability of accumulation introduced by local avalanche inputs at the glacier scale is not represented, and their impact on projections of future glacier mass changes has never been estimated. Here, we couple our gravitational snow redistribution model to the Open Global Glacier Model (OGGM), to explicitly account for avalanches. We then perform two sets of simulations: the ‘Control’ simulations, which correspond to the OGGM simulations without explicit representation of avalanches, and the ‘Avalanches’ simulations, where avalanches are represented using our gravitational snow redistribution model. For the ‘Avalanches’ simulations, we use the distributed avalanche correction factors calculated over 01/2000-12/2019 to adjust the initial precipitation along the glacier flowline before recalibrating the model parameters against the geodetic mass balance data from ref. 2 (Methods). We note that we did not update the avalanche correction factors over time (Methods). For the recalibration we use an informed three-step calibration procedure introduced in OGGM v1.6, starting with the recalibration of the precipitation correction factor. Given that avalanches tend to add mass to glaciers at the regional scale (Fig. 2), accounting for avalanches generally leads to lower precipitation correction values, with mean regional differences up to 0.1 for Central Europe (Supplementary Fig. 5). This is especially the case in regions with steep relief, which are more prone to avalanching (Fig. 2).

The glacier-wide mass balance is preserved given that we calibrate against the 01/2000-12/2019 geodetic mass balance from ref. 2. However, accounting for avalanches increases the spatial variability in the modelled mass balance of individual glaciers, which may also affect the mass balance in avalanche-free elevation bands due to the recalibration procedure (Fig. 3). These changes can produce different volume change trajectories between the ‘Control’ and ‘Avalanches’ simulations (Fig. 3). Argentière Glacier in the French Alps accumulates a considerable amount of mass by avalanching at the base of steep north-facing headwalls (Fig. 3a), and our modelled avalanche contribution of +17% agrees well with these remote sensing-based estimates. Our model shows that snow is redistributed from the upper headwalls onto the glacier surface in the accumulation area, and further down in the ablation area (Fig. 3a), with the recalibration leading to slightly more negative mass balance on the main glacier tongue (Fig. 3b). Despite the substantial contribution of avalanches to the accumulation there are limited differences in the future glacier volume change between the simulations with and without explicit representation of avalanches, which is consistent with the higher order ice flow simulations conducted by ref. 21.

Explicitly accounting for avalanches can have more effect on other glaciers and considerably change both their mass balance altitudinal patterns and their future volume changes. An example is Moore Glacier in New Zealand, located on a steep south-west facing slope (Fig. 3d–f). Most of the snow falling both in the accumulation and in the lower ablation zone of this glacier is redistributed to lower elevations including outside of the glacier, leading to 41% of the snow falling on this glacier to be removed by avalanches. This translates into a close-to-balance mass budget at the very top of the glacier, and a positive mass balance at the current glacier median-elevation (Fig. 3e, f). Accounting for this effect results in faster glacier retreat, with the glacier maintaining 57% less volume in the ‘Avalanches’ simulations than in the ‘Control’ simulations by the end of the century for the low emission (SSP1-2.6) scenarios (Fig. 3f). For the Black and Bazhin Glaciers in Alaska and South Asia West, a high glacier-wide avalanche contribution of +49% and +38%, respectively, translates in completely different behaviours (Fig. 3g–l). Indeed, on the one hand Black Glacier retreats much faster and reaches a volume 87% lower than in the ‘Control’ simulations by 2100 for SSP1-2.6, while Bazhin Glacier maintains 27% more volume. These differences in glacier projections for a similar glacier-wide avalanche contribution come from the local topography which can result in some cases in snow being redistributed from the accumulation zone into the ablation zone, thereby enhancing mass loss, and in other cases in the expansion of the accumulation zone by increasing the mass supply from higher elevation (Fig. 1). This also implies that while small avalanche-fed glaciers may be temporarily sheltered from atmospheric warming, they could also become very sensitive to it if rising temperatures affect the avalanche supply.

At the regional scale, explicitly accounting for avalanches reduces mass accumulation at high elevations, while the snow redistributed to low elevations lead to less negative mass balances (Fig. 4; Supplementary Fig. 7). The area-weighted mean mass balance is reduced by 0.05 m w.eq. yr−1 in the accumulation zone of Central European glaciers, and by 0.28 m w.eq. yr−1 in New Zealand (Fig. 4). Lower accumulation values are also visible in the other strongly avalanche-affected regions, such as the Caucasus and Middle East (0.11 m w.eq. yr−1 reduction), South Asia East (0.25) or the Low Latitudes (0.14), while the other regions are less affected (Supplementary Fig. 7). This reduction is most noticeable at the highest elevations, which also have a smaller area, so the effect on the glacier-wide mass balance is limited. The generally lower mass balance at high elevations in the strongly avalanche-affected regions highlight the sensitivity of the steep glaciers and glacierized headwalls to rising air temperatures. Indeed, these high-elevation glaciers have a low mass turnover because most of the deposited snow in these zones is gravitationally redistributed to lower elevations (Fig. 3). This enhanced sensitivity to warming at high elevations, highlighted by our regional-scale simulations, has been observed locally - for example in the Mont-Blanc massif, where retreating ice aprons are linked to increasing instability of steep mountain headwalls. At the glacier scale, reduced mass balance estimates observed when accounting for gravitational snow redistribution therefore indicate locations where glaciers are likely to shrink rapidly at high elevations, resulting in potential headwall instabilities due to the rapid exposure of permafrost (Fig. 3).

Longer subsistence of small glaciers with avalanches
Explicitly representing avalanches in OGGM has little effect on the regional-scale simulations of glacier volume and area changes from 2000 to 2100 (Fig. 5; Supplementary Figs. 8–25), with differences not exceeding 2%, even in the most avalanche-affected regions (Supplementary Fig. 8–25). This result arises because the mass balance model with or without avalanches is calibrated against the same geodetic mass balance data, and that the largest glaciers - generally less avalanche-affected (Fig. 2; Supplementary Fig. 1) - contribute the most to the regional volume change projections. Similar observations have been made in previous efforts to include sub-debris melt or calving in global glacier models.

However, given the influence of the avalanches on the altitudinal mass balance patterns, explicitly accounting for them can have a considerable effect at the scale of individual glaciers (Figs. 3, 5). When examining glaciers smaller than 1 km2 in the different RGI regions most affected by avalanches (Fig. 5; Supplementary Figs. 7–25), we find that accounting for avalanches leads to a slower volume loss than previously anticipated. Notably, the difference between the ‘Control’ and ‘Avalanches’ simulations peaks between 2020 and 2050, depending on the region and the SSP scenario. Furthermore, while volume losses under SSP3-7.0 and SSP5-8.5 are nearly identical by the end of the century with and without accounting for avalanches, in some regions small glaciers reach a different equilibrium under SSP1-2.6 with more volume preserved than when not accounting for avalanches.

In Central Europe, explicitly accounting for avalanches under SSP1-2.6 leads glaciers smaller than 1 km2 to retain 0.98 km3 of ice - three times more than in the ‘Control’ simulations (Fig. 5). For glaciers smaller than 0.5 km3, four times more ice is retained. In South Asia West, glaciers smaller than 1 km2 (0.5 km2) retain 19% (23%) more ice in the ‘Avalanches’ simulations (Supplementary Fig. 20). These results confirm the hypothesis of local-scale studies predicting the rapid and complete disappearance of small low-elevation glaciers, except for those strongly fed by avalanches. Here, we provide a first estimate of this effect across all mountain regions, which is important to consider given the role of small glaciers in downstream hydrology of mountain catchments.

Changing influence of avalanches as glaciers retreat
As glaciers retreat in a warming climate, the avalanche contribution is expected to be influenced by deglaciation of avalanche release and deposit slopes, glacier surface lowering, and changes in snowfall amounts including through the rise of the rain-snow transition line (Methods). The avalanche contribution may increase or decrease depending on mountain topography and future climate projections (Fig. 6, Supplementary Figs. 26–41). In Central Europe, glaciers are anticipated to retreat to the base of steep headwalls under all SSP scenarios, increasing the avalanche contribution from 11% in 2020 to 13–15% by 2100 (Fig. 6b). In New Zealand, the avalanche contribution depends strongly on the SSP scenario: it is expected to remain at its current level under SSP1-2.6, but decrease in higher emission scenarios, from 15% in 2020 to 3% in 2100 under SSP5-8.5 as the glacierized area shrinks (Fig. 6d–f). This likely reflects a dependence on the rising snow-line which, in high-warming scenarios, would rise above most of the relatively low summits in this range. In the Low Latitudes, glaciers will retreat to steeper, high-elevation slopes, causing relatively more accumulation to be removed by avalanches. As a result, the avalanche contribution is projected to decrease from −4% in 2020 to −10% in 2100 under SSP1-2.6 and from −4% to −30% under SSP5-8.5 (Fig. 6g–i). Changes in relative accumulation (in m w.eq. yr−1) appear to have a limited influence on the avalanche contribution, as these values remain consistent over time and across SSP scenarios in the different regions (Fig. 6c, f, i).

These simulations of avalanche contribution strongly depend on the topography of the glaciers and their surrounding catchments, and on the evolution of the glacier extents as simulated by our model. OGGM is a flowline model that represents the glacier evolution by elevation band (Methods) and requires an entire elevation band to be ice-free before removing it from the glacier extent, which can sometimes lead to inaccurate representations of the glacier geometry. Furthermore, while we are confident that our model correctly depicts the general regional trends, we note that in these simulations the avalanche contribution at the end of the century is dependent on the accurate geometric representation of the few glaciers remaining by 2100 in the regions most affected by glacier retreat. As a result, the uncertainty of these simulations is likely higher at the end of the century than the currently depicted range, which does not include uncertainties related to the glacier geometry. Despite these limitations, the contrasting regional patterns highlight the need to account for avalanches on individual glaciers, as overall glacier shrinkage will likely increase the role of these processes in the mass balance of glaciers across most mountain ranges.

Methods
Avalanche model description
Our avalanche model is based on the SnowSlide model, which has been widely used in glacio-hydrological studies of individual catchments. This model assumes an exponential decay of the maximum snow depth as a function of slope:
SNDmax = C e−αS (1)
Where SNDmax is the maximum snow depth at any given point of the digital elevation model (DEM) grid, S is the slope of the snow DEM (DEM with the snow layer on top), and C and α are the two parameters used to define this maximum snow depth, along with the minimum slope of 25° over which redistribution is allowed to occur. Once the maximum snow depth is defined, the excess snow is routed recursively using the multiple flow direction routing algorithm until no snow is redistributed anymore.

The large majority of studies that have applied this model have done it at hourly time-steps and using the original parameters proposed by ref. 31, originally obtained based on the comparison between end-of-season snow extents from Landsat ETM+ satellite images over the Waltzmann massif and the outputs of the SnowSlide model run at hourly time-steps and 30m spatial resolution. More recently, used the same maximum snow depth parametrization than 31 at 25, 50 and 100m resolution, but only using the minimum slope parameter of 25° for triggering, and reducing by 30% the maximum snow depth for pixels receiving snow, to mimic avalanche dynamics. They validated this approach against snow depth maps from repeat LiDAR acquisitions in the eastern Swiss Alps, which indicated a good agreement with the modelled snow depths.

For this study, we use values of 0.14 and 145 for the parameters α and C, respectively, along with the same parameterization of avalanche inertia introduced by ref. 32. We did not use a minimum slope threshold as we did not see any difference when including it. The values of α and C differ slightly from those originally reported by ref. 31 – 0.12 and 120, respectively. They are obtained based on a qualitative assessment of SnowSlide outputs in the Mt Blanc massif, and Argentière Glacier specifically, where we could directly compare our model outputs with distributed mass balance data from an ice flux inversion. This parametrization was further evaluated against Sentinel-1 avalanche deposits (Methods-Avalanche model evaluation).

We run the snow redistribution model at monthly time-steps and at a spatial resolution Δx (in m), varying as a function of glacier size S (in km2):
Δx = min(d1√S + d2, 200) (2)
Where the parameters d1 and d2 are set to 14 and 10, respectively. If the chosen spatial resolution is larger than 200m, we clip it to this value, so the spatial resolution effectively ranges between 10 and 200m. We tested the sensitivity of our results to the spatial resolution below.

Avalanche correction factors
We run SnowSlide independently for every month of the year using the 01/2000-12/2019 mean-monthly precipitation from W5E5v2.0 data as inputs, along with a density of 200 kg m−3 to convert the solid precipitation inputs to a snow height. SnowSlide is run over the whole glacier and its surroundings within a buffer of 80 grid cells (with a spatial resolution varying from 10m to 200m) around the glaciers, which we assume to encompass the entire glacier catchment. We then sum the monthly snow maps and take the ratio of snow accumulation before and after redistribution to compute a distributed avalanche correction factor kava. The values within the glacier outlines from the Randolph Glacier Inventory (RGI) version 6.0 are then aggregated:
1. Over the entire glacier to obtain a glacier-wide avalanche correction factor kava,gl:
2. Along the glacier flowline to obtain an avalanche correction factor for each elevation band kava,i following 60.

We note that the use of RGI 6.0 outlines, for example due to potential misalignments, can lead to uncertainties in glacier future projections and also most likely in estimates of the avalanche contribution. At the regional scale the uncertainties are small relative to other contributions such as from the choice of climate models, but for individual glaciers these estimates should be interpreted with caution.

Recalibration of the glacier model
The avalanche correction factors kava,i vary along the glacier flowline and can be applied for each glacier as multiplicative factors to the monthly solid precipitation PSolid(z) (in kg m−2 month−1) to compute the mean monthly mass balance Bi(z) using the temperature index model from OGGM building from 27, 29, 63 and expressed as:
Bi(z) = kava,i(z) × kp × PSolid(z) − df × max(ti + T(z) − T Melt, 0) (3)
Where Ti(z) is the mean monthly air temperature (in °C) at altitude z using a −6.5°C km−1 lapse rate and T Melt is the monthly mean air temperature above which ice melt is assumed to occur (−1°C per default). In OGGM precipitation is directly taken from the nearest climate reanalysis cell, without any correction for elevation, contrary to the temperature. The precipitation is considered to fall in solid (liquid) form when the mean monthly temperature is below 0°C (above 2°C), with a linear transition in-between.

The precipitation correction factor kp (-), the degree-day factor df (in kg m−2 K−1 month−1) and the temperature bias ti (in K) are three parameters that are calibrated for each glacier against the 01/2000-12/2019 average geodetic mass balance from ref. 2 and using W5E5v2.0 data. This calibration uses the informed three-steps approach introduced in OGGM v1.6, which is run until the geodetic mass balance is matched following these steps:
1. the three parameters are prescribed some allowed ranges (kp ∈ [0.1, 10], df ∈ [0.33, 33] kg m−2 K−1 day−1, ti ∈ [−8, 8] K),
2. kp is obtained from the glacier’s winter precipitation from W5E5v2.0 based on a logarithmic relation derived by ref. 29,
3. df is fixed at 5 kg m−2 K−1 day−1, ti is calibrated to match the geodetic mass balance of each glacier, and its median value per W5E5v2.0 grid point is attributed to each glacier in the grid point,
4. kp is recalibrated allowing it to vary between +/- 20% of its original value, then df and then ti are recalibrated.

When accounting for avalanching, we directly use the W5E5v2.0 monthly solid precipitation multiplied by the avalanche correction factor (kava,i × PSolid(z)) and the pre-calibrated kp, df and ti. We then recalibrated these parameters starting from step 4. In practice, introducing avalanches has two main interconnected effects after calibration: first, accumulation patterns along the glacier are modified according to the avalanche model; second, the initial precipitation correction factor, originally estimated without accounting for avalanche contributions, is recomputed based on this new information. In the cases where no initial precipitation factor can be found within the specified +/− 20% bounds, df and ti are changing as well, as described in step 4.

Following these different steps, we obtain two mass balance models, both calibrated against the 01/2000-12/2019 geodetic mass balance for all glaciers in the world: the ‘Control’ model which does not explicitly account for avalanches, and the ‘Avalanches’ model, where avalanche contribution is explicitly accounted for using the kava,i avalanche correction factors.

Avalanche contribution to glacier accumulation
We compute the net glacier-wide avalanche contribution to accumulation Avgl, expressed as a percentage:
Avgl = (SNDwithava,gl − SNDnoava,gl) / SNDwithava,gl × 100 (4)
Where SNDnoava,gl and SNDwithava,gl are the snow depths before and after avalanche redistribution, integrated over the entire glacier scale.
This net glacier-wide avalanche contribution can be split between the positive avalanche contribution where snow is added to the glacier (Av+gl), and the negative avalanche contribution where snow is removed from the glacier (Av−gl):
Avgl = Av+gl + Av−gl (5)
Similarly, we compute the net region-wide avalanche contribution for each of the 19 RGI regions and for the entire globe by summing the snow accumulation before and after avalanches for all glaciers of the regions:
Avreg = (∑glaciers(SNDwithava,gl − SNDnoava,gl)) / ∑glaciers SNDwithava,gl × 100 (6)

Comparison to glacier metrics
We compare the net, positive and negative glacier-wide mean avalanche contributions calculated over the period 01/2000-12/2019 with a series of glacier metrics. The average slope along the flowline, size and mean elevation are derived from the RGI 6.0 glacier attributes, the geodetic mass balance is obtained from the dataset by ref. 2 and the percentage of debris cover is obtained from the global dataset by ref. 64.

Glacier projections
We simulate the evolution of all glaciers in the world for the period 2000–2100 with both the ‘Control’ and the ‘Avalanches’ models. As climate data for 2020–2100 we use all the available models out of the 18 climate models of the Coupled Model Intercomparison Project (CMIP) Phase 6 for the three different Shared Socioeconomic Pathways (SSP1-2.6, SSP3-7.0 and SSP5-8.5). These climate data are downscaled and bias-corrected to the W5E5 data following the procedure described in ref. 65. We retain all the glaciers for which at least 75% of the SSP/GCM combinations ran without errors, which corresponds to 98.6% of the RGI6.0 glaciers (97.7% of the total glacier area; Supplementary Table 2). For these projections we take the avalanche correction factor calculated over the period 01/2000-12/2019 and do not update it with time, because of the high associated computational costs. We did however test the influence of updating this avalanche correction factor every 20 years for 1255 glaciers with a glacier-wide avalanche correction factor higher than 1.2 or lower than 0.8 (corresponding to an avalanche contribution higher than 17% or lower than -25%), randomly selected in 15 of the RGI regions for the ipsl-cm6a-lr_r1i1p1f1 GCM under SSP1-2.6 (Supplementary Figs. 42–43). We found the differences between the avalanche scenarios with or without temporal update of the avalanche correction factor to be in most cases negligible relative to the differences between the avalanche and the control scenarios. This suggests that although glacier surface lowering (creating steeper and taller headwalls) and a rising rain-snow transition may affect avalanche-driven mass redistribution in the future, these influences are likely secondary compared to the main climatic controls on glacier evolution. In the few cases where the temporal update of avalanches made a noticeable difference like for the glaciers with a positive avalanche contribution in Alaska or North Asia (Supplementary Fig. 43), it increased the gap between the control and avalanche scenarios. This suggests that not updating the avalanche contribution, as done in this study, gives a lower-bound estimate of its long-term influence on glacier evolution.

Avalanche model evaluation
While SnowSlide has been extensively used and validated, and is considered to give a realistic estimate of avalanche redistribution for hydrological models, we extend these validation efforts to glacierized regions using remote sensing observations.
We validate our avalanche model against avalanche deposit extents derived from Sentinel-1 images over a 5-year period for the Mt Blanc, Everest and Hispar regions. While the data derived from Sentinel-1 is limited to deposit extents and does not give any information on redistributed volume or mass, we show that the deposits from our model occur at the same elevation range as the deposits derived from Sentinel-1 images, albeit covering a generally smaller area, especially in the Hispar region (Supplementary Figs. 49–51). These results also highlight a strong correlation between the observed and modelled avalanche deposit areas for individual glaciers (R2 ≥ 0.73), which confirms that our model adequately represents the spatial domain of the observed avalanche activity.

Direct measurements of avalanche contribution to glacier accumulation would be more useful but it is notoriously difficult to measure from remote sensing or field data. For this, we mostly rely on the validation efforts by ref. 32, although they were conducted for non-glacierized regions. Some studies have estimated the local contribution of avalanches to glacier mass balance using Ground Penetrating Radar measurements, but highlighted a strong interannual and spatial variability, making it challenging to attribute the variability to avalanches and to obtain a glacier-wide contribution estimate. We compared the modelled glacier-wide precipitation correction factor with the value obtained by ref. 21 from an inversion of the ice flux constrained by remote sensing observations for Argentière Glacier, finding a good agreement between the two with the measured and modelled avalanche correction factor kava,gl being equal to 1.2 and 1.3, respectively, for the glacier without its North tributary, and 1.1 and 1.3 with all the tributaries.

Ultimately, we compared the modelled and measured annual mass balance for all measurements available over the period 01/2000-12/2019 in the global database from the World Glacier Monitoring Service. 9646 annual point mass balance measurements were retained, acquired across 90 glaciers distributed across all RGI regions except Arctic Canada South, the Russian Arctic, North Asia, the Caucasus and Middle East, South Asia West and the Subantarctic and Antarctic Islands. Including avalanches in the mass balance model had a negligible impact on the comparison with the mass balance measurements, with a RMSE of 1.15 m w.eq. yr−1 for the Control scenario and 1.17 m w.eq. yr−1 for the Avalanches scenario, and the same R2 (Supplementary Fig. 52). This is not surprising given that due to representativity and safety considerations, the huge majority of mass balance measurements are conducted away from any avalanche deposits and the long-term monitoring sites have a limited estimated avalanche contribution, with 90% of these glaciers having an avalanche contribution between −2% and +10% (Supplementary Fig. 53).

Sensitivity of the model to DEM spatial resolution
Given the computational constraints of simulating snow redistribution on all glaciers of the world and surrounding catchments, we use a spatial resolution that varies with glacier size, between 10 and 200m. Coarser spatial resolutions translate into less steep slopes, and ultimately in less snow being redistributed by our gravitational snow redistribution model. In the mountainous region of South-East Switzerland, showed that SnowSlide simulations at 25m resolution resulted in the best agreement with LiDAR measurements of snow volumes, but they still found similar snow redistribution patterns when running their simulations at 50 or 100m, with significant improvement from the simulations without gravitational snow redistribution.

Here, we test the sensitivity of the avalanche contribution and glacier volume changes to the spatial resolution of the DEMs by running our model at 50m and 100m resolution for Central Europe (Supplementary Table 3; Supplementary Figs. 44–46). We find the regional avalanche contribution to be similar in all three cases, both for the historical period and for the future simulations, with the results from the reference simulations being more similar to those from the 50m resolution simulations (Supplementary Figs. 44–46). Such findings are also true at the scale of small individual glaciers (Supplementary Figs. 47–48). Furthermore, a lower resolution lead to lower avalanche contribution and to a less strong effect of avalanching on the future projections of small glaciers. While this does not influence the findings of this study, it is worth noting that our variable spatial resolution may lead to an underestimation of the avalanche contribution to larger glaciers relative to the very small ones.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "This is a valuable contribution that explores the role of avalanche to glacier mass gain (or loss) at the global scale – the first time such an analysis has been applied to glaciers globally. The authors couple an open source (OGGM) glacier model to a gravitational snow redistribution model (Snow Slide) to explore glacier mass balance with and without modelled avalanche input. Model parameters were calibrated to match the geodetic mass balance developed by Hugonnet et al (2021) and results are also compared to glaciological mass balance data. Future projections were then made for using three SSPs exploring how avalanche accumulation may affect glacier evolution in a warming world. This is a high-quality contribution, it is well written, well cited, and contains a number of quality Figures to support findings. I enjoyed reading this contribution and have only minor comments below.\nSpecific comments to each section follow:\nTitle: Perhaps the term ‘variable' rather than ‘significant would be better in the title? While the percentage contribution of avalanche input is indeed significant in some glacier regions, many other regions had contributions of less than 1%. When considering Supplementary Fig 6, one might determine that across the global sample the avalanche contribution was not ‘significant’ but quite significant particularly in alpine settings.\nLine 34: warming - minor typo ‘warmin’\nFig 1: This is a nice Figure, but the blue lines in 1a are a bit confusing. Does the blue line indicate a hypothetical glacier area?\nLine 66: Reference to the ‘temperature bias’ can this be explained more clearly. Is this referring to a known bias associated with the gridded temperature data, to the known variability of temperature with elevation, or, be a reference to previous reports of elevation dependant warming? A few additional sentences to explain exactly what temperature bias authors are referring would be appreciated.\nLine 91: Model simulations have been set-up to align with calendar years. Was there any consideration to using the hydrological year – even though northern and southern hemisphere will be opposing? Starting simulations in January means that the first ablation (accumulation) season is half over in the southern (northern) hemisphere respectively, although equatorial glaciers will be different again. Perhaps this is why the start of the time series look a bit odd in Fig 6 for the NZ region? But maybe this does not matter depending on model spin-up processes?\nLine 108 (plus Sup): It is acknowledged that the Randolf Glacier Inventory is the ‘standard’ for global glacier outlines. However, there exists some rather inaccurate outlines for some glaciers. For example, the RGI outline of Rolleston Glacier (NZ) is at best a poorly aligned representation of the glacier back in the 1980-90s. It is significantly oversize (and misaligned). Since at least 2010 the entire portion of the RGI outline to the south-east is no longer connected to the current glacier, it is well separated by bedrock and not considered part of the glacier. The slopes above this now non-existent part of the glacier are prone to avalanche deposition, so one might expect this will create a positive bias in results. Without in-depth knowledge of all glacier regions, one wonders how many other small glaciers suffer the same poor alignment of RGI outlines, which will consequently impact the accuracy of the geodetic mass balance as produced by Hugonnet et al and results of this study. Did the authors undertake spot-checks of the glacier areas? I appreciate that this is perhaps beyond the scope of the analysis, but such errors are potentially perpetuating from one global study to the next.\nLines 136-138: The results relating to where on the glacier avalanche deposition occurs is very interesting and ideally warrants further investigation/discussion. Avalanche deposition onto the ablation zone creates a variable range of feedbacks, especially in relation to albedo/melting. Avalanche deposition onto an ablation area may increase or decrease albedo depending on surface type (ice/debris) and avalanche deposition is known to be of higher density, thus more mass. Could the authors say a little more about the differences between glaciers where avalanche deposition is primarily occurring in the accumulation area of glacier compared to avalanche deposition the occurs in glacier ablation areas?\nLine 168: It would seem that a number of regions required almost no recalibration of the precipitation factor, which would indicate that at these sites the role of avalanche input (either positive or negative is negligible). Results shown in Sup Fig 5 which are close to null, align with results in Fig 2 for areas where avalanche processes are negligible. Once might conclude that in general it is glaciers in mountainous terrain for which avalanche deposition is important. Perhaps this finding could be stated more explicitly.\nLine 175-200: The individual glacier case studies are interesting and a great way to tease out specific examples. Could the authors include a Figure (possibly in supplementary) that supports Fig 3 so readers can quickly ascertain exactly where in the regional mountain ranges the various case studies are located.\nLine 305-307: As glaciers thin slope angle tends to increase in the headwall region this may be an additional ‘future’ feedback.\nLine 311-315: The Southern Alps of NZ are relatively low elevation compared to mountain ranges around the globe, perhaps this lack of ‘real-estate’ at higher elevation contributes to the finding that avalanche contribution strongly depends on the SSP. As temperatures warm the snow/rain threshold will likely be off the top of many mountains in NZ.\nLines 350-364: The authors have highlighted some interesting results; this is a great study. However, I think that considering results only in terms of ‘geopolitical’ regions limits the findings. What are the patterns for cirque glaciers globally compared to valley glaciers (regardless of region)? What about debris-covered glaciers? It appears that the Himalayan glaciers create a bias, I can think of many debris-covered glacier for which avalanche accumulation is not an important contribution. How do glaciers that have significant avalanche deposition in their ablation areas compare to those where avalanche is confirmed to the upper accumulation areas? These are all very interesting spatial questions for which the authors have results, but which are not discussed?\nMethods section: This was well written and clear to follow. I found that questions I had on reading the core paper e.g. impact of DEM resolution, role of precipitation factors etc were all answered in due course in this section.\nSupplementary: The supplementary is substantial, but all the material is well presented and strongly supports the main document. It was good to see all the regional graphs, although seeing some Figures/results focused on glacier characteristics rather than glacier regions would have been a great addition.", "reviewer_2": "Review of: \"Substantial contribution of avalanches to glacier accumulation in the 21st century\", M. Kneib et al.\nGeneral comment\nThis paper explores the contribution of gravitational snow redistribution to glacier accumulation, through an extensive study of the world glaciers by coupling the simple avalanche model SnowSlide to the global glacier mass balance model Open Global Glacier Model (OGGM). The evolution of this contribution throughout the 21st century is presented, using different emission scenarios. The authors quantify the global avalanche contribution to glacier accumulation (3%) and erosion (1%), with stronger contributions in some regions like New Zealand and for smaller glaciers. The spatial variability of glacier mass balance is increased by avalanches, with local impacts on individual glacier evolutions. At regional scale, avalanches have little impact on the evolution of glacier volumes in the 21st century, but have a significant impact for small glacier with longer subsistence. The avalanche contribution itself has varying evolutions in the 21st century depending on regions.\nThis study brings novel and useful insights, by explicitly adding the contribution of snow avalanches to glacier mass balance modelling. This explicit modelling improves existing empirical parameterizations which can’t account for the non-linearity of this process. The worldwide application of the method throughout the 21st century provides very interesting results for the cryosphere community. The methods used are generally sound, properly justified and clearly presented. The paper is well written and clearly exposes its main conclusions. I only have a few minor comments and questions regarding the methods in snow redistribution, detailed hereafter. Note that my review may be incomplete since my expertise does not cover glacier modelling itself.\nSpecific minor comments\n• Parameterization of SnowSlide\nThe calibration of the snow holding depth following the slope differs from the original one. Why? Can you indicate how this function was calibrated, e.g. where and using what data? Do you think the calibration could be different depending on the regions, their topography and their snow climatology?\n• Influence of varying spatial resolution\nAs you mentioned, smaller grid spacing on smaller glaciers leads to more avalanching than on larger glaciers. The study of sensitivity to spatial resolution is useful and interesting. However, it could be supplemented by showing sensitivity for the most extreme cases (the smallest glaciers), since these are an important focus of the main results. How would the results change for a very small glacier (e.g. Moore Glacier in Fig. 3) where spatial resolution is around 20 m, if this resolution was the same as for the largest glaciers (200 m)? Would the outcomes for these very small glaciers change significantly?\nTechnical comments and typos\n• The readability of some figures, particularly Fig. 2, can be improved: the text is sometimes hard to read at 100% zoom.\n• The following comment is more of a note, as I assume it is of limited significance in the context of global modelling and long-term climate projections.\nAspect dependency: I understand that no snowpack model is used prior to SnowSlide, which is, of course, justified in the context of the study. But since the snow depth used by SnowSlide is a proxy for snowfall, doesn't it lead to avalanching in sun-exposed faces being relatively overestimated versus shaded faces, where less snow melt and compaction could occur before avalanching?\n• L. 34: warming\n• L. 64-65: I assume you mean: snow redistribution by wind\n• L. 129: “absolute” or relative?\n• L. 272: Fig. 7-24", "reviewer_3": "I read this model study with great interest and was impressed by the authors’ efforts to quantify the avalanche contribution for glaciers worldwide. This work represents a valuable contribution, and the dataset provided offers an excellent basis for further analyses of the factors driving avalanche contributions to glaciers and how these may evolve under a changing climate. This study constitutes a very interesting first step. However, I have several comments and concerns that I believe should be addressed before publication.\nMain comments:\nAvalanche process representation: Avalanche contributions to glaciers are often linked to enhanced snow accumulation resulting from preferential deposition on surrounding slopes. As far as I understand, this process is not explicitly accounted for in the study. Instead, the avalanche contribution appears to be modeled primarily as a function of topography (slope angle and elevation of adjacent terrain and the glacier terrain), which may not fully capture the underlying processes. This simplification could lead to an underestimation of avalanche effects, particularly for many small glaciers. That said, I recognize the necessity of using a simplified approach when modeling avalanche contributions at a global scale, and the study represents an important step forward. However, I would encourage the authors to include a critical discussion of this limitation and its potential impact on the model results.\nModel time step: Coming from a background in small-scale studies, I wonder why the authors chose a relatively coarse model time step of one month. How does this affect the representation of precipitation distribution and the partitioning between solid and liquid precipitation?\nSection order: I was surprised by the ordering of the sections. I would have preferred reading the methods before the results, as this structure provides the necessary foundation for understanding what the authors have actually done.\nMinor comments:\nIntroduction: Please provide an overview of the range of avalanche contributions estimated in previous studies for individual glaciers.\nResults – topographic relationships: It would be valuable to analyze the relationship between avalanche contribution (in percentage) and the aspect, slope, and elevation difference of surrounding slopes and the glacier slopes that contribute to the overall mass to glaciers. This could help assess which glaciers may benefit more from avalanches in the future, separated by these topographic characteristics.\nClarity in results and discussion: I appreciate the analysis of potential future effects of avalanches on glacier mass balance. However, the presentation of the results and discussion is currently quite wordy and sometimes confusing. I recommend rewriting this section to present the information more clearly, emphasizing the interdependencies between processes, climate change, and glacier size evolution.\nPrecipitation processes: Please discuss the effects of missing precipitation processes and their strong spatial variability, particularly in high-elevation and steep mountain terrain. Include how these omissions could influence your projections for the future.\nMethods section:\nVariable grid size: The use of a variable grid size depending on glacier size is understandable, but I suggest discussing the potential effects on model results. Changing the grid size affects the representation of topography, slope angles, and therefore snow redistribution, which impact the model outputs and introduces a bias for smaller versus larger glaciers.\nMinimum slope angle of 25°: Please provide reasons for not using the minimum slope angle of 25°.\nFigures:\n• The fonts in many of the figures are very small and difficult to read; please increase the font size for better clarity. Additionally, the quality of Figure 1 appears somewhat low and could be improved.\n• Figure 1: I do not clearly understand the distinction between “1. Snow removal” and “4. Snow redistribution.” In the figure, “Snow removal” also seems to indicate redistribution from very steep, high-elevation slopes to flatter, low-elevation areas of the glacier. Please revise this figure to make the processes and their differences more clear."}, "author_response": "We’d like to thank reviewer 1 for their comments and thorough review, this is really appreciated.\nSpecific comments to each section follow:\nTitle: Perhaps the term ‘variable' rather than ‘significant would be better in the title? While the percentage contribution of avalanche input is indeed significant in some glacier regions, many other regions had contributions of less than 1%. When considering Supplementary Fig 6, one might determine that across the global sample the avalanche contribution was not ‘significant’ but quite significant particularly in alpine settings.\nWe’ve changed the title to ‘Topographically-controlled contribution of avalanches to glacier mass balance in the 21st century’.\nLine 34: warming - minor typo ‘warmin’\nGood catch. Modified (L36).\nFig 1: This is a nice Figure, but the blue lines in 1a are a bit confusing. Does the blue line indicate a hypothetical glacier area?\nWe have modified this figure to distinguish the avalanche release, transition and deposit areas and added a legend for the glacier outlines and mountain ridges. We have also changed the color of the glacier outlines & crevasses so that it is not confused with the surface mass balance ‘Control’ line:\nFigure 1: Impact of avalanches occurring at the surface of a hypothetical mountain glacier. a\nHypothetical glacier with avalanche release, transition and deposit areas in light blue, purple and red,\nrespectively. b Corresponding effect on the annual mass balance patterns with elevation. Avalanches\noccurring within the glacier boundaries (1-3 & 6) lead to additional accumulation on the avalanche\ndeposits (3) while the removal of snow in the release areas locally reduces the surface mass balance\n(1). Avalanches from outside (4-5) the glacier boundaries, from seasonally snow-covered slopes or\nother glaciers, lead to a local increase in mass balance.\nLine 66: Reference to the ‘temperature bias’ can this be explained more clearly. Is this referring to a\nknown bias associated with the gridded temperature data, to the known variability of temperature\nwith elevation, or, be a reference to previous reports of elevation dependant warming? A few\nadditional sentences to explain exactly what temperature bias authors are referring would be\nappreciated.\nThis ‘temperature bias’ is a parameter of the mass balance model used in OGGM - an additive\nparameter to the temperature supposed to correct potential offsets in the climate reanalysis\ntemperature series, indeed for the reasons the reviewer mentions (the bias parameter is temporally\nconstant however, i.e. it cannot account for elevation dependant warming). We’ve now clearly\nspecified that this is a model parameter (L70-72):\n‘This precipitation correction factor, together with a temperature bias parameter and the melt factor\nused in global-scale temperature-index models, needs to be calibrated on a glacier-by-glacier basis28.’\nWe note that this is explained in detail in the Methods (Recalibration of the glacier model). We\ntherefore do not insist on this here in the Introduction.\nLine 91: Model simulations have been set-up to align with calendar years. Was there any\nconsideration to using the hydrological year – even though northern and southern hemisphere will\nbe opposing? Starting simulations in January means that the first ablation (accumulation) season is\nhalf over in the southern (northern) hemisphere respectively, although equatorial glaciers will be\ndifferent again. Perhaps this is why the start of the time series look a bit odd in Fig 6 for the NZ\nregion? But maybe this does not matter depending on model spin-up processes?\nThanks for this thoughtful input. The period 01/2000-12/2019 had to be the same as the glacier model\ncalibration period, which directly has to be identical to the geodetic mass balance estimates from\nHugonnet et al. (2021), as is done for other global modeling studies (Rouce et al., 2023; Zekollari et\nal., 2025). This was necessary to stay consistent between the period over which the avalanche inputs\nwere calculated, and the glacier model calibration period. As a result, that indeed means that we start\nin January both with the glacier and avalanche simulations but we would expect this effect to be minor\nand just affect the very first few months - note that the avalanche contribution was then averaged\nover this whole period, so these first few months are unlikely to have any influence on the model runs.\nNew Zealand does show an increase in glacier area at the start of the historical period (which we guess\nis what was meant when pointing at the oddness in Fig. 6), contrary to the other regions (Fig. 6; Fig.\nS25-S40). This could be due to an initial shock due to the absence of model spinup, but this is unlikely\nas this is not visible for the other regions. Perhaps this could be related to the small glacier advance\nthat occurred in the early 2000s due to some cooler temperatures at the time (as suggested by\nMackintosh et al., 2017). This also is visible in both satellite remote sensing (Zemp et al., 2025) and in\nsitu mass balance observations (Sirguey et al., 2016; Cullen et al., 2019).\nWe have now specified in the text that 01/2000-01/2019 is the model calibration period to clarify this\nchoice (L95-97):\n‘We use a simple gravitational snow redistribution model31,32 run with W5E5v2.0 mean monthly\nclimate inputs33 over the glacier model calibration period 01/2000-12/2019 to quantify the\ncontribution of avalanches to the annual glacier accumulation (Methods).’\nLine 108 (plus Sup): It is acknowledged that the Randolf Glacier Inventory is the ‘standard’ for global\nglacier outlines. However, there exists some rather inaccurate outlines for some glaciers. For\nexample, the RGI outline of Rolleston Glacier (NZ) is at best a poorly aligned representation of the\nglacier back in the 1980-90s. It is significantly oversize (and misaligned). Since at least 2010 the\nentire portion of the RGI outline to the south-east is no longer connected to the current glacier, it is\nwell separated by bedrock and not considered part of the glacier. The slopes above this now non-\nexistent part of the glacier are prone to avalanche deposition, so one might expect this will create a\npositive bias in results. Without in-depth knowledge of all glacier regions, one wonders how many\nother small glaciers suffer the same poor alignment of RGI outlines, which will consequently impact\nthe accuracy of the geodetic mass balance as produced by Hugonnet et al and results of this study.\nDid the authors undertake spot-checks of the glacier areas? I appreciate that this is perhaps beyond\nthe scope of the analysis, but such errors are potentially perpetuating from one global study to the\nnext.\nWe fully acknowledge the limitations of RGI 6.0, which we depended on in this study given that the\ngeodetic mass balance estimates needed for the calibration were produced specifically for this\ninventory (Hugonnet et al., 2021). Rolleston is a typical example of these limitations and this was one\nof the reasons that we did not attempt to validate our model at this site, despite it being one of the\nrare glaciers with estimates of avalanche contribution (Purdie et al., 2015). As a result, the glacier-\nscale simulations need to be interpreted with caution and we showed them only to indicate the\nvariability in avalanche influence between different glaciers (Fig. 3). Currently all global-scale glacier\nmodeling studies rely on the RGI (until now RGI 6.0, given that the main glacier products - ice thickness,\ngeodetic mass balance - have not yet been adapted to RGI 7.0).\nWe did not conduct spot-checks ourselves as this was indeed far beyond the scope of this global-scale\nstudy, and such limitations from the RGI 6.0 and their influence on model results have been described\nalready in the scientific literature (e.g. Hartl et al., 2025). However, we refer the reviewer to the work\nby Aguayo et al. (2024) who assessed the uncertainty of OGGM in the Chilean Andes resulting from\nthe choice of glacier inventory. They compared the future glacier projections obtained with RGI 6.0\nand the with the Chilean national inventory, but also with various ice thickness estimates or climate\nscenarios and demonstrated that ‘While the glacier inventory and ice thickness source showed overall\ndifferences close to 10 %, the different climate alternatives showed differences of more than 50 % for\nsolid precipitation, for example.’ (Aguayo et al., 2024).\nUltimately, we note that this study aims to provide a first estimate of the influence of avalanches on\nglacier mass balance at the global scale. It is only a first step that will undoubtedly be far from\nperfect but has the advantage of putting some numbers on a process that was until now completely\nignored.\nWe however agree that these limitations need to be clearly stated and discussed and have now\nexpanded their description in the Methods (L488-492):\n‘We note that the use of RGI 6.0 outlines, for example due to potential misalignments, can lead to\nuncertainties in glacier future projections and also most likely in estimates of the avalanche\ncontribution61. At the regional scale the uncertainties are small relative to other contributions such as\nfrom the choice of climate models62, but for individual glaciers these estimates should be interpreted\nwith caution.’\nLines 136-138: The results relating to where on the glacier avalanche deposition occurs is very\ninteresting and ideally warrants further investigation/discussion. Avalanche deposition onto the\nablation zone creates a variable range of feedbacks, especially in relation to albedo/melting.\nAvalanche deposition onto an ablation area may increase or decrease albedo depending on surface\ntype (ice/debris) and avalanche deposition is known to be of higher density, thus more mass.\nCould the authors say a little more about the differences between glaciers where avalanche deposition is\nprimarily occurring in the accumulation area of glacier compared to avalanche deposition the occurs in glacier ablation areas?\nThank you. Please see our answer to the comment below for the comparison of glaciers with\navalanche contribution in the ablation versus accumulation area.\nFor the feedbacks related to albedo and density, this is an excellent point, which however cannot be\naddressed with the simple model that we used: the OGGM mass balance model that we use does\nnot include any representation of snow metamorphism, nor therefore albedo or density - see\nMethods. This choice was made to be able to get a first global estimate of the avalanche\ncontribution and its primary effect on glacier evolution, but considerably simplifies the related\nprocesses. Almost all global glacier models tend to rely on models with a limited number of\nparameters to reduce the equifinality when calibrating against the geodetic mass balance only\n(Schuster et al., 2023), and we therefore stayed with this choice to be able to compare the runs with\navalanches to the OGGM standard simulations.\nWe are currently working on a field-based study where we aim to disentangle these different\nprocesses occurring on one avalanche deposit of Argentière Glacier (French Alps), but this is a whole\ndifferent study conducted at a completely different spatial scale. For this paper, we have however\nnow acknowledged these limitations in the main text (L407-414):\n‘A major challenge, however, lies in calibrating and validating gravitational snow redistribution\nmodels in glacierized areas due to the scarcity of remote sensing and field-based observations.\nExpanding such datasets will be crucial to refine future model estimates (Methods), particularly to\ncapture feedbacks related to avalanche deposition, such as effects on albedo and snow density, and\nto better understand how avalanche activity depends on the characteristics of release areas (aspect,\nslope, or hanging glaciers).’\nLine 168: It would seem that a number of regions required almost no recalibration of the\nprecipitation factor, which would indicate that at these sites the role of avalanche input (either\npositive or negative) is negligible. Results shown in Sup Fig 5 which are close to null, align with\nresults in Fig 2 for areas where avalanche processes are negligible. Once might conclude that in\ngeneral it is glaciers in mountainous terrain for which avalanche deposition is important. Perhaps\nthis finding could be stated more explicitly.\nAbsolutely. We have insisted on this in the main text (L175-178):\n‘Given that avalanches tend to add mass to glaciers at the regional scale (Fig. 2), accounting for\navalanches generally leads to lower precipitation correction values, with mean regional differences\nup to 0.1 for Central Europe (Supplementary Figure 5). This is especially the case in regions with steep\nrelief, which are more prone to avalanching (Fig. 2).’\nLine 175-200: The individual glacier case studies are interesting and a great way to tease out specific\nexamples. Could the authors include a Figure (possibly in supplementary) that supports Fig 3 so\nreaders can quickly ascertain exactly where in the regional mountain ranges the various case studies\nare located.\nWe have added a Supplementary Figure 6 showing the precise location of these glaciers:\nSupplementary Figure 6: Location of the four focus glaciers. The glaciers are indicated by a red\ntriangle and the blue outlines in the background show the RGI 6.0 glacier outlines in the corresponding\nregions: a Argentière Glacier (RGI6.0-11.03638; 45.951°N, 6.985°E), b Moore Glacier (RGI6.0-\n18.02066; -43.898°N, 169.786°E), c Black Glacier (RGI6.0-01.13640; 61.085°N, -143.831°E) and d\nBazhin Glacier (RGI6.0-14.20029; 35.226°N, 74.637°E).\nLine 305-307: As glaciers thin slope angle tends to increase in the headwall region this may be an\nadditional ‘future’ feedback.\nAbsolutely, this is partly what we implied when mentioning ‘glacier surface lowering’. We have\nspecified this in more details in the Methods (L575-579):\n‘This suggests that although glacier surface lowering (creating steeper and taller headwalls) and a\nrising rain-snow transition may affect avalanche-driven mass redistribution in the future20,49, these\ninfluences are likely secondary compared to the main climatic controls on glacier evolution. ’\nLine 311-315: The Southern Alps of NZ are relatively low elevation compared to mountain ranges\naround the globe, perhaps this lack of ‘real-estate’ at higher elevation contributes to the finding that\navalanche contribution strongly depends on the SSP. As temperatures warm the snow/rain threshold\nwill likely be off the top of many mountains in NZ.\nThis is an excellent suggestion, that we have added to the main text (L329-334):\n‘In New Zealand, the avalanche contribution depends strongly on the SSP scenario: it is expected to\nremain at its current level under SSP1-2.6, but decrease in higher emission scenarios, from 15% in\n2020 to 3% in 2100 under SSP5-8.5 as the glacierized area shrinks (Fig. 6d-f). This likely reflects a\ndependence on the rising snow-line which, in high-warming scenarios, would rise above most of the\nrelatively low summits in this range.’\nLines 350-364: The authors have highlighted some interesting results; this is a great study. However,\nI think that considering results only in terms of ‘geopolitical’ regions limits the findings. What are the\npatterns for cirque glaciers globally compared to valley glaciers (regardless of region)? What about\ndebris-covered glaciers? It appears that the Himalayan glaciers create a bias, I can think of many\ndebris-covered glacier for which avalanche accumulation is not an important contribution. How do\nglaciers that have significant avalanche deposition in their ablation areas compare to those where\navalanche is confirmed to the upper accumulation areas? These are all very interesting spatial\nquestions for which the authors have results, but which are not discussed?\nWe really appreciate these inputs, that refer to the potential of this new dataset to explore\nrelationships and test the influence of avalanches for various types of glaciers. There would be a lot\nto do, but at this stage we mainly wanted to provide a first quantification geared towards regional\nvalues, along with open-access to the code and dataset, to allow for more in-depth analysis at\nvarious spatial scales. We therefore refrained from adding too much of these analyses to the text.\nThe different points that were raised are however definitely worth exploring and we provide short\nanswers to these questions below. This as led us to make a few modifications of the text that we also\nlist here:\n● Debris-covered glaciers: many regions show a significant positive correlation between debris\ncover and avalanching (Fig. S2-S4), but there is a lot of spread so this obviously doesn’t\nexplain the complete signal. This is true not just for the HMA regions, but also for Alaska,\nWestern Canada and USA, Central Europe or New Zealand, especially when considering only\nthe glaciers receiving a net accumulation surplus by avalanching (Fig. S3). We totally agree\nthat many debris-covered glaciers are not avalanche-fed and vice versa. We show however\nthat there is a certain tendency of debris-covered glaciers to also be avalanche-fed (which\ndoes not indicate causality!). These correlations actually likely come from the fact that\ndebris-covered glaciers and avalanche-fed glaciers alike need to be surrounded by steep\nterrain for debris OR snow to be deposited on their surface. But we are indeed looking at\ntwo different kinds of processes here. We have highlighted this in the description of these\ncorrelations (L141-143):\n‘This confirms past observations conducted at the scale of High Mountain Asia41, which likely\nindicate that debris-covered and avalanche-fed glaciers exist in a similar topographic setting,\nwith steep surrounding slopes.’\n● Cirque versus valley glaciers: This would be an excellent comparison, albeit a not so\nstraightforward one at the global scale, which was the focus of this study. Indeed, the\nclassification does not exist for all 200,000+ RGI glaciers. It would be possible only for a\nsubset of glaciers for which these attributes have been indicated in the GLIMS database. As\nto whether these glaciers would be influenced differently by avalanches definitely calls for\nmore analysis which is however beyond the scope of this study.\n● Avalanches affecting the ablation versus accumulation areas: yes, this is also a good point,\nwhich could potentially explain the different behaviours of Black and Bazhin Glaciers despite\ntheir similar avalanche contribution, in Figure 3. This is also what we hinted at in the main\ntext (L211-216):\n‘These differences in glacier projections for a similar glacier-wide avalanche contribution\ncome from the local topography which can result in some cases in snow being redistributed\nfrom the accumulation zone into the ablation zone, thereby enhancing mass loss, and in\nother cases in the expansion of the accumulation zone by increasing the mass supply from\nhigher elevation (Fig. 1). ’\nHowever looking at this at the large scale is far from trivial for three main reasons: 1/ many\nglaciers receive avalanches both in their ablation and accumulation areas, some also lose\nsnow by avalanching from both. A classification - while definitely possible - is thus not\nstraightforward and would require defining a specific metric to characterize this ratio. 2/ The\naccumulation-ablation area ratio is changing with time, sometimes rapidly with climate\nchange, and different mountain regions are responding more or less fast to this climate\nchange. This makes it challenging to define a reference time at which to make this\ncomparison. 3/ Accounting for avalanches leads to a change of definition of ‘accumulation’\nand ‘ablation’ zone. Removing snow from the accumulation area may lead to parts of this\narea to have a negative mass balance while an avalanche deposit in what we traditionally\ncategorize as ‘ablation’ area could have positive mass balance (Fig. 3).\nDespite the underlying challenges, we do believe that such analysis would be worthwhile and should\nbe undertaken, but we felt that this was beyond the scope of this study. We have however now\nacknowledged this potential for new analysis in the main text (L388-392):\n‘This study also paves the way for more detailed sub-regional studies, for example examining links\nbetween avalanches and debris cover (Supplementary Figures 2-3), differences in avalanche effects\nacross glacier types (valley, cirque), or the impacts of avalanches occurring in accumulation versus\nablation zones (Fig. 3).’\nMethods section: This was well written and clear to follow. I found that questions I had on reading\nthe core paper e.g. impact of DEM resolution, role of precipitation factors etc were all answered in\ndue course in this section.\nThank you.\nSupplementary: The supplementary is substantial, but all the material is well presented and strongly\nsupports the main document. It was good to see all the regional graphs, although seeing some\nFigures/results focused on glacier characteristics rather than glacier regions would have been a great\naddition.\nThanks, it’s actually really exciting to explore this new dataset and look at individual glaciers.\nHowever, as explained in the comments above, we decided to stick to the regional estimates, which\nare also less uncertain. Figure 3 was the exception to show a few examples of where the regional\npatterns could come from. This paper will be complemented by a more process-focused study to be\nsubmitted in the coming months.\nReviewer #1 (Remarks on code availability):\nI fully support the authors providing the code and this is likely something I will visit and explore in\ndue course, and potentially utilise for post-graduate projects. However, I have not (as yet) used the\nOGGM, so did not check the code.\nThanks, we’ll be very happy if this code can be used by others! The documentation is not (yet) to the\nstandard of OGGM’s but hopefully it is not too difficult to handle.\nReviewer #2 (Remarks to the Author):\nReview of: \"Substantial contribution of avalanches to glacier accumulation in the 21st century\", M.\nKneib et al.\nGeneral comment\nThis paper explores the contribution of gravitational snow redistribution to glacier accumulation,\nthrough an extensive study of the world glaciers by coupling the simple avalanche model SnowSlide\nto the global glacier mass balance model Open Global Glacier Model (OGGM). The evolution of this\ncontribution throughout the 21st century is presented, using different emission scenarios. The\nauthors quantify the global avalanche contribution to glacier accumulation (3%) and erosion (1%),\nwith stronger contributions in some regions like New Zealand and for smaller glaciers. The spatial\nvariability of glacier mass balance is increased by avalanches, with local impacts on individual glacier\nevolutions. At regional scale, avalanches have little impact on the evolution of glacier volumes in the\n21st century, but have a significant impact for small glacier with longer subsistence. The avalanche\ncontribution itself has varying evolutions in the 21st century depending on regions.\nThis study brings novel and useful insights, by explicitly adding the contribution of snow avalanches\nto glacier mass balance modelling. This explicit modelling improves existing empirical\nparameterizations which can’t account for the non-linearity of this process. The worldwide\napplication of the method throughout the 21st century provides very interesting results for the\ncryosphere community. The methods used are generally sound, properly justified and clearly\npresented. The paper is well written and clearly exposes its main conclusions. I only have a few\nminor comments and questions regarding the methods in snow redistribution, detailed hereafter.\nNote that my review may be incomplete since my expertise does not cover glacier modelling itself.\nWe’d like to thank reviewer 2 for their comments and thorough review, this is really appreciated.\nSpecific minor comments\n• Parameterization of SnowSlide\nThe calibration of the snow holding depth following the slope differs from the original one. Why?\nCan you indicate how this function was calibrated, e.g. where and using what data? Do you think the\ncalibration could be different depending on the regions, their topography and their snow\nclimatology?\nThanks for this comment. Indeed, the parameters we used (0.14; 145) vary slightly from the original\nones (0.12; 120) - or at least these are the values that we found when digitizing the Fig. 2b plot from\nBernhardt & Schulz (2010), originally obtained based on the comparison between end-of-season snow\nextents from Landsat ETM+ satellite images over the Watzmann massif.\nFirst, we note that the differences between the maximum snow depth as a function of slope for the\ndifferent parameter pairs ([0.14; 145] VS [0.12; 120]) are very small, as indicated in the figure below\n(Fig. R1), which translates to the avalanche deposit extents and volumes being also similar for both\nparmetrizations.\nFigure R1: Maximum snow depth as a function of slope relationships using two different sets of\nparameters for the snow redistribution model.\nThe parametrization that we use tends to remove a bit more snow from slopes in the [10; 40]° range\nthan the original parametrization. This parametrization was obtained based on a qualitative\nassessment of SnowSlide outputs in the Mt Blanc massif, and Argentière Glacier specifically, where we\ncould directly compare the OGGM-SnowSlide outputs with distributed mass balance data from an ice\nflux inversion 21. There, we found a good agreement with the runs at monthly time scale and the\nadditional parametrization of inertia from Quéno et al. (2024). In our study we further evaluated the\navalanche extents obtained with this parametrization against the Sentinel-1 avalanche deposit extents\nfrom Kneib et al. (2024a), as described in the methods.\nWe have now mentioned this in the methods (L458-466):\n‘For this study, we use values of 0.14 and 145 for the parameters a and C, respectively, along with the\nsame parameterization of avalanche inertia introduced by 32. We did not use a minimum slope\nthreshold as we did not see any difference when including it. The values of a and C differ slightly from\nthose originally reported by 31 – 0.12 and 120, respectively. They are obtained based on a qualitative\nassessment of SnowSlide outputs in the Mt Blanc massif, and Argentière Glacier specifically, where we\ncould directly compare our model outputs with distributed mass balance data from an ice flux inversion\n21. This parametrization was further evaluated against Sentinel-1 avalanche deposits (Methods -\nAvalanche model evaluation).’\nBased on this, we could also envision that the SnowSlide parameters may need to be recalibrated on\na regional basis, most likely more because of the snow climatology than the topography. This would\nbe a very interesting test, but which would require the appropriate data for these calibration runs,\nideally along the lines of what Quéno et al. (2024) did using repeat snow depth measurements over\nlarge areas using data from an airborne LiDAR. This is now mentioned in the discussion (L407-414):\n‘A major challenge, however, lies in calibrating and validating gravitational snow redistribution\nmodels in glacierized areas due to the scarcity of remote sensing and field-based observations. Expanding such\ndatasets will be crucial to refine future model estimates (Methods), particularly to capture feedbacks\nrelated to avalanche deposition, such as effects on albedo and snow density, and to better understand\nhow avalanche activity depends on the characteristics of release areas (aspect, slope, or hanging\nglaciers).’\n• Influence of varying spatial resolution\nAs you mentioned, smaller grid spacing on smaller glaciers leads to more avalanching than on larger\nglaciers. The study of sensitivity to spatial resolution is useful and interesting. However, it could be\nsupplemented by showing sensitivity for the most extreme cases (the smallest glaciers), since these\nare an important focus of the main results. How would the results change for a very small glacier\n(e.g. Moore Glacier in Fig. 3) where spatial resolution is around 20 m, if this resolution was the same\nas for the largest glaciers (200 m)? Would the outcomes for these very small glaciers change\nsignificantly?\nThanks for this comment. Actually Supplementary Figures 44 & 45 d-f show the effect of lowering the\nspatial resolution to 50 and 100 m for all the Central Europe glaciers smaller than 1 km2 (so with a\nresolution of 24 m or higher). These results show significant changes especially when going from 50 m\nto 100 m resolution, less when going from the original resolution to 50 m (Fig. 5), but the overall\npatterns are maintained.\nWe get similar results when looking at this on a glacier-by-glacier basis, as shown below for two small\nglaciers in Central Europe, where the altitudinal mass balance and future volume changes are\nparticularly affected when going from 50 m to 100 m resolution, less so in the step from the original\nresolution to a resolution of 50 m. We have added this as two new figures in the supplementary\ninformation to complement the regional-scale observations, and added a sentence on this topic to the\nmain text:\nMain text (L639-646):\n‘Here, we test the sensitivity of the avalanche contribution and glacier volume changes to the spatial\nresolution of the DEMs by running our model at 50 m and 100 m resolution for Central Europe\n(Supplementary Table 3; Supplementary Figures 44-46). We find the regional avalanche contribution\nto be similar in all three cases, both for the historical period and for the future simulations, with the\nresults from the reference simulations being more similar to those from the 50 m resolution simulations\n(Supplementary Figures 44 & 46). Such findings are also true at the scale of small individual glaciers\n(Supplementary Figures 47-48). ’\nSupplementary figures:\nSupplementary Figure 47: Influence of DEM spatial resolution on volume change of Tré-les-Eaux\nGlacier. a Tré-les-Eaux Glacier in Central Europe with altitudinal mass balance and 21st century volume\nchanges b-c at the original 20 m resolution, d-e at 50 m resolution and f-g at 100 m resolution. The left\npanel shows a 3D view of the glacier with the average avalanche correction factor (k ) over the period\nava\n01/2000-12/2019 indicated by the red-white-blue shade wrapped over the glacier extents from the\nRandolph Glacier Inventory (RGI) 6.039. The background topography is from the DEM used to run the\nOGGM simulations. The numbers indicate the glacier-wide avalanche contribution and the glacier\narea. The central panels show the altitudinal mass balance for the ‘Control’ and ‘Avalanches’\nsimulations, with the red dashed line corresponding to a mass balance value of 0 m w.eq. yr-1. The right\npanels show the past and future volume changes of these glaciers as modelled using OGGM. All\npercentages are given relative to the initial glacier volume in 2000 in the Control scenario. The black\nline corresponds to the historical period over which the mass balance model was calibrated using\nW5E5v2.0 data33. The colored lines show the median future projections for different SSP scenarios, and\nthe shaded areas indicate the 25th-75th percentile range. The different curves were smoothed using a\n5-year rolling mean.\nSupplementary Figure 48: Influence of DEM spatial resolution on volume change of Arête Midi-Plan.\na Arête Midi-Plan in Central Europe with altitudinal mass balance and 21st century volume changes b-\nc at the original 20 m resolution, d-e at 50 m resolution and f-g at 100 m resolution. The left panel\nshows a 3D view of the glacier with the average avalanche correction factor (k ) over the period\nava\n01/2000-12/2019 indicated by the red-white-blue shade wrapped over the glacier extents from the\nRandolph Glacier Inventory (RGI) 6.039. The background topography is from the DEM used to run the\nOGGM simulations. The numbers indicate the glacier-wide avalanche contribution and the glacier\narea. The central panels show the altitudinal mass balance for the ‘Control’ and ‘Avalanches’\nsimulations, with the red dashed line corresponding to a mass balance value of 0 m w.eq. yr-1. The right\npanels show the past and future volume changes of these glaciers as modelled using OGGM. All\npercentages are given relative to the initial glacier volume in 2000 in the Control scenario. The black\nline corresponds to the historical period over which the mass balance model was calibrated using\nW5E5v2.0 data33. The colored lines show the median future projections for different SSP scenarios, and\nthe shaded areas indicate the 25th-75th percentile range. The different curves were smoothed using a\n5-year rolling mean.\nTechnical comments and typos\n• The readability of some figures, particularly Fig. 2, can be improved: the text is sometimes hard to\nread at 100% zoom.\nWe have now increased the font sizes of figures 5, 6 and figure 2:\nWe have also updated figures 1 & 3 for readability, increasing the font sizes and adding a divergent\ncolorbar to figure 3:\n• The following comment is more of a note, as I assume it is of limited significance in the context of\nglobal modelling and long-term climate projections.\nAspect dependency: I understand that no snowpack model is used prior to SnowSlide, which is, of\ncourse, justified in the context of the study. But since the snow depth used by SnowSlide is a proxy\nfor snowfall, doesn't it lead to avalanching in sun-exposed faces being relatively overestimated\nversus shaded faces, where less snow melt and compaction could occur before avalanching?\nThis is an excellent point which we totally agree with. As the reviewer states, there was not really\nany alternative to SnowSlide to make these calculations at the global scale but we agree that this\nwould be an important point to consider at smaller spatial scales. We have now added a sentence\nrelated to this point in the last section of the main text (L407-414):\n‘A major challenge, however, lies in calibrating and validating gravitational snow redistribution\nmodels in glacierized areas due to the scarcity of remote sensing and field-based observations.\nExpanding such datasets will be crucial to refine future model estimates (Methods), particularly to\ncapture feedbacks related to avalanche deposition, such as effects on albedo and snow density, and\nto better understand how avalanche activity depends on the characteristics of release areas (aspect,\nslope, or hanging glaciers).’\n• L. 34: warming\nGood catch. Modified (L36).\n• L. 64-65: I assume you mean: snow redistribution by wind\nIndeed! Added (L70).\n• L. 129: “absolute” or relative?\nGood catch, we’ve only left ‘the avalanche contribution’ as the definition speaks for itself (L136).\n• L. 272: Fig. 7-24\nThanks, this has been adapted to the new figure indexing.\nReviewer #3 (Remarks to the Author):\nI read this model study with great interest and was impressed by the authors’ efforts to quantify the\navalanche contribution for glaciers worldwide. This work represents a valuable contribution, and the\ndataset provided offers an excellent basis for further analyses of the factors driving avalanche\ncontributions to glaciers and how these may evolve under a changing climate. This study constitutes\na very interesting first step. However, I have several comments and concerns that I believe should be\naddressed before publication.\nMany thanks for this thorough review and comments which we’ve tried to address as best as we\ncould.\nMain comments:\nAvalanche process representation: Avalanche contributions to glaciers are often linked to enhanced\nsnow accumulation resulting from preferential deposition on surrounding slopes. As far as I\nunderstand, this process is not explicitly accounted for in the study. Instead, the avalanche\ncontribution appears to be modeled primarily as a function of topography (slope angle and elevation\nof adjacent terrain and the glacier terrain), which may not fully capture the underlying processes.\nThis simplification could lead to an underestimation of avalanche effects, particularly for many small\nglaciers. That said, I recognize the necessity of using a simplified approach when modeling avalanche\ncontributions at a global scale, and the study represents an important step forward. However, I\nwould encourage the authors to include a critical discussion of this limitation and its potential\nimpact on the model results.\nThanks for this comment. Indeed our approach relies solely on precipitation and topography to\nestimate the avalanche contribution of each individual glacier. We do not account for any wind\nredistribution, which could locally enhance or reduce the effective snow depth, and therefore affect\nthe avalanche supply. Wind speed and direction being a very localized process, it is challenging to\ntake into account beyond the glacier/catchment scale, and such a modelling endeavour would likely\nbe very difficult to calibrate/validate (Voordendag et al., 2024).\nThis was suggested in the original manuscript (L103-107):\n‘We note that this model does not represent wind-driven snow redistribution, which is a particularly\ndifficult process to constrain given its dependence on an accurate representation of wind fields on a\nglacier, and despite recent advances in quantifying and modelling this effect at relatively small\nspatial scales32,36–38.’\nWe have also now specifically mentioned this point in the discussion (L420-432):\n‘Our results underscore the need for an improved representation of local processes affecting the\nspatial variability of mass balance in glacier models, which, beyond avalanches, include frontal\nablation of marine- and lake-terminating glaciers, variable melt under debris, variable precipitation\nwith elevation or snow redistribution by wind. Further research is needed to understand how these\nprocesses interact30,32,38,53 and how they can be effectively integrated into global glacier models6,8,30.\nFor example, wind-driven snow redistribution may locally enhance or reduce the effective snow depth\non slopes on or in the vicinity of glaciers, and therefore affect the avalanche supply. These process\nrepresentations, combined with the recent developments in remote sensing data availability and\nassimilation by global glacier models and the anticipated capability of running these simulations in a\ndistributed way, will pave the way for a new generation of models capable of accurately projecting\nthe future of every individual glacier across Earth’s mountain ranges.’\nModel time step: Coming from a background in small-scale studies, I wonder why the authors chose\na relatively coarse model time step of one month. How does this affect the representation of\nprecipitation distribution and the partitioning between solid and liquid precipitation?\nThis is a very good point. This choice was driven mainly by computational constraints and also because\n1 month is the time-step at which OGGM is run, which is also the temporal resolution of the vast\nmajority of current global glacier models (Zekollari et al., 2025; Wimberly et al., 2025). Therefore, we\nalso used the same precipitation partitioning scheme than OGGM, which is applied to monthly\nprecipitation & temperature values.\nSchuster et al. (2023) specifically tested running OGGM at daily resolution versus monthly and did not\nfind a systematic effect on glacier volume changes, as both models were calibrated against the same\ndata. In that study, more than the precipitation partitioning, the authors emphasized the resulting\nchanges in monthly melt as being responsible for the increase in the spread of volume changes - see\na full discussion of this effect in paragraph 4.1.2 of Schuster et al. (2023).\nBased on these considerations, we decided to stay consistent with the reference OGGM runs and keep\na one-month time step for our own simulations. A higher temporal resolution for SnowSlide would\nanyway not have been realistic at this spatial scale.\nSection order: I was surprised by the ordering of the sections. I would have preferred reading the\nmethods before the results, as this structure provides the necessary foundation for understanding\nwhat the authors have actually done.\nThis order of the sections is imposed by the Journal, and we can therefore not change it.\nMinor comments:\nIntroduction: Please provide an overview of the range of avalanche contributions estimated in\nprevious studies for individual glaciers.\nWe have now added this (L50-53):\n‘In the few locations where on-glacier avalanching has been estimated, it is one of the main sources\nof mass inputs, with estimated avalanche contributions ranging between 9% and 92% of the annual\naccumulation13,14,19–22. ’\nResults – topographic relationships: It would be valuable to analyze the relationship between\navalanche contribution (in percentage) and the aspect, slope, and elevation difference of\nsurrounding slopes and the glacier slopes that contribute to the overall mass to glaciers. This could\nhelp assess which glaciers may benefit more from avalanches in the future, separated by these\ntopographic characteristics.\nThis comment hints towards attempting to derive an avalanche contribution from topographic\ncharacteristics only, and without any mass redistribution modeling. We consider this to be beyond the\nscope of this study and give a detailed summary of our reasons below.\nWhile we agree that this is a worthwhile objective, we object to doing this with our results for the\nfollowing reasons:\n1. These are only model results, and not measurements or estimates of the actual avalanche\ncontribution. While such a comparison would be extremely worthwhile for measured\navalanche contributions, we expect that if we were to conduct it with our modelling results,\nit would be biased with the SnowSlide assumptions (homogeneous snow depth with aspect,\ndensity assumptions, spatial resolution used…)\n2. Conducting such an analysis would be equivalent to trying to simplify the already quite\nsimplistic avalanche estimation scheme. The SnowSlide scheme already accounts for slopes\nand precipitation, and can already be run relatively efficiently at the large scale as\ndemonstrated by this study. Our recommendation would therefore be to directly use this\nredistribution scheme to estimate the avalanche contribution to glacier mass balance rather\nthan relying on a proxy.\n3. We expect the contributions from the different topographic factors to be difficult to\ndisentangle, and that there will not be a single metric that can explain the diversity of\navalanche contributions. This is based on a preliminary analysis of the avalanche contribution\nas a function of the R-index, which is the proportion of slopes steeper than 30° within a\nglacierized catchment. This metric has been used as a proxy for the potential of avalanche\ncontribution in past studies, and proven to be one of the best predictors for the avalanche\ncontribution (Hughes, 2008; Kneib et al., 2024a). We compared the SnowSlide avalanche\ncontribution to this R index for 478 glaciers for which the R-index had been calculated in a\nrecent study by Kneib et al. (2024a). We find that while a certain value of R-index is a necessary\ncondition for a glacier to have a given avalanche contribution, it does not explain the whole\nvariability in avalanche contributions (Fig. R1). These results are similar to those of Kneib et\nal. (2024a) conducted with mapped Sentinel-1 avalanche deposits in the Mt. Blanc, Everest\nand Hispar regions.\nFigure R1: Avalanche contribution (from this study) as a function of R-index (from Kneib et al., 2024a)\nfor 478 glaciers in the Mt Blanc, Hispar and Everest regions. Each circle represents one glacier and\nthe circle size is scaled to the glacier size.\nBased on this and the fact that such an analysis would require a substantial amount of work to be\nconducted at the global scale, we have not taken this any further.\nClarity in results and discussion: I appreciate the analysis of potential future effects of avalanches on\nglacier mass balance. However, the presentation of the results and discussion is currently quite\nwordy and sometimes confusing. I recommend rewriting this section to present the information\nmore clearly, emphasizing the interdependencies between processes, climate change, and glacier\nsize evolution.\nWe have been through the entire Results & Discussion and simplified the text as much as we could.\nPlease refer to the tracked-changes version of the manuscript to see all the minor changes that were\nmade.\nPrecipitation processes: Please discuss the effects of missing precipitation processes and their strong\nspatial variability, particularly in high-elevation and steep mountain terrain. Include how these\nomissions could influence your projections for the future.\nThis is of course particularly relevant, not only for avalanches but for all the accumulation processes.\nSince all we did in this study was to add the representation of a process in the OGGM mass balance\nmodel, we do not expect the uncertainty of the precipitation to change between the model with and\nwithout explicit representation of avalanches. Regional glacier models correct precipitation based on\nglacier-specific calibration (Schuster et al., 2023), and this correction is a major source of uncertainty\nin regional simulations and is the main cause of inter-model runoff differences (Marzeion et al.,\n2020; Wimberly et al., 2025).\nWe have now specifically mentioned the lack of representation of variable precipitation (L420-432):\n‘Our results underscore the need for an improved representation of local processes affecting the\nspatial variability of mass balance in glacier models, which, beyond avalanches, include frontal\nablation of marine- and lake-terminating glaciers, variable melt under debris, variable precipitation\nwith elevation or snow redistribution by wind. Further research is needed to understand how these\nprocesses interact30,32,38,53 and how they can be effectively integrated into global glacier models6,8,30.\nFor example, wind-driven snow redistribution may locally enhance or reduce the effective snow depth\non slopes on or in the vicinity of glaciers, and therefore affect the avalanche supply. These process\nrepresentations, combined with the recent developments in remote sensing data availability and\nassimilation by global glacier models and the anticipated capability of running these simulations in a\ndistributed way, will pave the way for a new generation of models capable of accurately projecting\nthe future of every individual glacier across Earth’s mountain ranges. ’\nMethods section:\nVariable grid size: The use of a variable grid size depending on glacier size is understandable, but I\nsuggest discussing the potential effects on model results. Changing the grid size affects the\nrepresentation of topography, slope angles, and therefore snow redistribution, which impact the\nmodel outputs and introduces a bias for smaller versus larger glaciers.\nPlease see our detailed response to the comment from Reviewer 2.\nMinimum slope angle of 25°: Please provide reasons for not using the minimum slope angle of 25°.\nWhen testing the model in various mountain ranges we did not see any difference when including\nthis minimum slope threshold, which is why we did not use this condition. We have now specified\nthis in the text (L459-460):\n‘We did not use a minimum slope threshold as we did not see any difference when including it.‘\nFigures:\n• The fonts in many of the figures are very small and difficult to read; please increase the font size\nfor better clarity. Additionally, the quality of Figure 1 appears somewhat low and could be improved.\nWe have now increased the font size of figures 2, 3, 5 and 6. For figure 1, please see our detailed\nresponse to the comment below.\n• Figure 1: I do not clearly understand the distinction between “1. Snow removal” and “4. Snow\nredistribution.” In the figure, “Snow removal” also seems to indicate redistribution from very steep,\nhigh-elevation slopes to flatter, low-elevation areas of the glacier. Please revise this figure to make\nthe processes and their differences more clear.\nWe have modified this figure to distinguish the avalanche release, transition and deposit areas and\nadded a legend for the glacier outlines and mountain ridges. We have removed the process\ndescriptions “1. Snow removal” and “Snow redistribution” and used a more explicit terminology to\ndescribe the elements of the figure:\nFigure 1: Impact of avalanches occurring at the surface of a hypothetical mountain glacier. a\nHypothetical glacier with avalanche release, transition and deposit areas in light blue, purple and red.\nb Corresponding effect on the annual mass balance patterns with elevation. Avalanches occurring\nwithin the glacier boundaries (1-3 & 6) lead to additional accumulation on the avalanche deposits (3)\nwhile the removal of snow in the release areas reduces the surface mass balance (1). Avalanches from\noutside (4-5) the glacier boundaries, from seasonally snow-covered slopes or other glaciers, lead to a\nlocal increase in mass balance.\nREFERENCES\nAguayo, R., Maussion, F., Schuster, L., Schaefer, M., Caro, A., Schmitt, P., Mackay, J., Ultee, L., Leon-\nMuñoz, J., & Aguayo, M. (2024). Unravelling the sources of uncertainty in glacier runoff projections\nin the Patagonian Andes (40–56° S). The Cryosphere, 18(11), 5383–5406. https://doi.org/10.5194/tc-\n18-5383-2024\nBernhardt, M., & Schulz, K. (2010). SnowSlide: A simple routine for calculating gravitational snow\ntransport. Geophysical Research Letters, 37(11). https://doi.org/10.1029/2010GL043086\nCullen, N. J., Gibson, P. B., Mölg, T.,Conway, J. P., Sirguey, P., & Kingston,D. G. (2019). The influence\nof weather systems in controlling mass balance in the Southern Alps of New Zealand. Journal of\nGeophysical Research:Atmospheres, 124, 4514–4529. https://doi.org/10.1029/2018JD030052\nHartl, L., Schmitt, P., Schuster, L., Helfricht, K., Abermann, J., & Maussion, F. (2025). Recent\nobservations and glacier modeling point towards near-complete glacier loss in western Austria\n(Ötztal and Stubai mountain range) if 1.5 °C is not met. The Cryosphere, 19(3), 1431–1452.\nhttps://doi.org/10.5194/tc-19-1431-2025\nHughes, Philip. D. (2008). Response of a montenegro glacier to extreme summer heatwaves in 2003\nand 2007. Geografiska Annaler: Series A, Physical Geography, 90(4), 259–267.\nhttps://doi.org/10.1111/j.1468-0459.2008.00344.x\nHugonnet, R., McNabb, R., Berthier, E., Menounos, B., Nuth, C., Girod, L., Farinotti, D., Huss, M.,\nDussaillant, I., Brun, F., & Kääb, A. (2021). Accelerated global glacier mass loss in the early twenty-\nfirst century. Nature 2021 592:7856, 592(7856), 726–731. https://doi.org/10.1038/s41586-021-\n03436-z\nJouberton, A., Shaw, T. E., Miles, E., Kneib, M., Fugger, S., Buri, P., McCarthy, M., Kayumov, A.,\nNavruzshoev, H., Halimov, A., Kabutov, K., Homidov, F., & Pellicciotti, F. (2025). Snowfall decrease in\nrecent years undermines glacier health and meltwater resources in the Northwestern Pamirs.\nCommunications Earth & Environment, 6(1), 691. https://doi.org/10.1038/s43247-025-02611-8\nKneib, M., Dehecq, A., Brun, F., Karbou, F., Charrier, L., Leinss, S., Wagnon, P., & Maussion, F. (2024).\nMapping and characterization of avalanches on mountain glaciers with Sentinel-1 satellite imagery.\nThe Cryosphere, 18(6), 2809–2830. https://doi.org/10.5194/tc-18-2809-2024\nMackintosh, A. N., Anderson, B. M., Lorrey, A. M., Renwick, J. A., Frei, P., & Dean, S. M. (2017).\nRegional cooling caused recent New Zealand glacier advances in a period of global warming. Nature\nCommunications, 8(1), 14202. https://doi.org/10.1038/ncomms14202\nMarzeion, B., Hock, R., Anderson, B., Bliss, A., Champollion, N., Fujita, K., Huss, M., Immerzeel, W.\nW., Kraaijenbrink, P., Malles, J. H., Maussion, F., Radić, V., Rounce, D. R., Sakai, A., Shannon, S., van\nde Wal, R., & Zekollari, H. (2020). Partitioning the Uncertainty of Ensemble Projections of Global\nGlacier Mass Change. Earth’s Future, 8(7), e2019EF001470. https://doi.org/10.1029/2019EF001470\nPurdie, H., Rack, W., Anderson, B., Kerr, T., Chinn, T., Owens, I., & Linton, M. (2015). The impact of\nextreme summer melt on net accumulation of an avalanche fed glacier, as determined by ground‐\npenetrating radar. Geografiska Annaler: Series A, Physical Geography, 97(4), 779–791.\nhttps://doi.org/10.1111/geoa.12117\nQuéno, L., Mott, R., Morin, P., Cluzet, B., Mazzotti, G., & Jonas, T. (2024). Snow redistribution in an\nintermediate-complexity snow hydrology modelling framework. The Cryosphere, 18(8), 3533–3557.\nhttps://doi.org/10.5194/tc-18-3533-2024\nSchuster, L., Rounce, D. R., & Maussion, F. (2023). Glacier projections sensitivity to temperature-\nindex model choices and calibration strategies. Annals of Glaciology, 1–16.\nhttps://doi.org/10.1017/aog.2023.57\nSirguey, P., Still, H., Cullen, N. J., Dumont, M., Arnaud, Y., and Conway, J. P.: Reconstructing the mass\nbalance of Brewster Glacier, New Zealand, using MODIS-derived glacier-wide albedo, The\nCryosphere, 10, 2465–2484, https://doi.org/10.5194/tc-10-2465-2016, 2016.\nWimberly, F., Ultee, L., Schuster, L., Huss, M., Rounce, D. R., Maussion, F., Coats, S., Mackay, J., &\nHolmgren, E. (2025). Inter-model differences in 21st century glacier runoff for the world’s major\nriver basins. The Cryosphere, 19(4), 1491–1511. https://doi.org/10.5194/tc-19-1491-2025\nZemp, M., Jakob, L., Dussaillant, I., Nussbaumer, S. U., Gourmelen, N., Dubber, S., A, G., Abdullahi, S.,\nAndreassen, L. M., Berthier, E., Bhattacharya, A., Blazquez, A., Boehm Vock, L. F., Bolch, T., Box, J.,\nBraun, M. H., Brun, F., Cicero, E., Colgan, W., … Zheng, W. (2025). Community estimate of global\nglacier mass changes from 2000 to 2023. Nature, 639(8054), 382–388.\nhttps://doi.org/10.1038/s41586-024-08545-z\nZekollari, H., Schuster, L., Maussion, F., Hock, R., Marzeion, B., Rounce, D. R., Compagno, L., Fujita,\nK., Huss, M., James, M., Kraaijenbrink, P. D. A., Lipscomb, W. H., Minallah, S., Oberrauch, M., van\nTricht, L., Champollion, N., Edwards, T., Farinotti, D., Immerzeel, W., … Sakai, A. (2025). Glacier\npreservation doubled by limiting warming to 1.5°C versus 2.7°C. Science, 388(6750), 979–983.\nhttps://doi.org/10.1126/science.adu4675\nWe would like to thank all three reviewers for their excellent, thorough and swift comments on our\nmanuscript. This has been a really nice peer-review process."}, {"round": 2, "reviewer_comments": {"reviewer_1": "Dear Marin,\nThank-you for your considered responses to all the reviewers feedback on your manuscript. The detail provided and edits made really strengthen this contribution. This is a great piece of work, which I am sure will generate much interest. I look forward to seeing this work published.\nKind regards\nHeather", "reviewer_2": "The authors provided thorough answers to the reviewers in general, and to my comments in particular. The additions and modifications improve the quality of the manuscript. I recommend publication. Congratulations for this very interesting study.", "reviewer_3": "All of my questions and comments have been successfully addressed, and the manuscript has been thoroughly and appropriately revised. I am satisfied with the changes made, and I recommend the manuscript for publication."}, "author_response": ""}]

Decision letter