We extracted both supra-glacial debris cover and clean-ice outlines from Scherler et al. (2018) for our glacier sample from 2014 to compare these results of our regional study with those from the global study. We found that a large portion of these glaciers in the Greater Caucasus are covered by supra-glacial debris cover. While for two regions the results match within the uncertainty (western Greater Caucasus: 22.4±6.4 % vs. 23.7 %; Elbrus: 4.6±6.6 % vs. 7.7 %), our values are lower than the results of Scherler et al. (2018) for the other two regions (10.1±6.2 % vs. 30.8 % for central, 49.2±5.7 % vs. 84.9 % for eastern; Fig. S3). These differences can mostly be explained by (i) the RGI v6, used by Scherler et al. (2018), which is characterized by some inconsistent co-registration for the Greater Caucasus region that stems from the use of inadequately orthorectified satellite imagery to generate the inventory in contrast to the improved orthorectification of the Landsat L1T data (Fig. S4a), and (ii) the RGI v6, containing nominal glaciers (i.e. ellipses around glacier label points) for the Greater Caucasus region, which originates from the use of the World Glacier Inventory (WGI; Haeberli et al., 1989) to fill gaps with no data for earlier versions of the RGI (Pfeffer et al., 2014). According to Scherler et al. (2018), all nominal glaciers were classified as debris-covered (Fig. S4b). We note that the scope of the study by Scherler et al. (2018) was an automatized global assessment of supra-glacial debris cover from optical satellite data, without correcting any outlines in the RGI.

Our results are in agreement with studies assessing changes in supra-glacial debris cover for smaller areas in the Greater Caucasus. For example, Stokes et al. (2007) calculated that supra-glacial debris cover increased by 3 %–6 % (∼0.20 % yr −1 ) between 1985 and 2000 on several glaciers in the central Greater Caucasus. On individual glaciers, supra-glacial debris cover increased by 25 %–30 % (e.g. Shkhelda) in the same period. We found an increase in supra-glacial debris cover, from 21.3±6.0 % to 30±5.8 % (∼0.65 % yr −1 ), for Shkhelda Glacier between 1986 and 2000. Popovnin et al. (2015), reported a supra-glacial debris cover increase from 2 % to 13 % (∼0.23 % yr −1 ) between 1968 and 2010, based on direct field monitoring for Djankuat Glacier (northern slope of the central Greater Caucasus). This compares with our result of an area increase in supra-glacial debris cover for Djankuat Glacier from 2.6±6.9 % to 9.8±6.1 % (∼0.25 % yr −1 ) between 1986 and 2014. This difference can be explained in that the detailed field measurement would have picked out smaller spots of debris cover which were beyond the resolution of the satellite imagery. Lambrecht et al. (2011) estimated that the supra-glacial debris cover distribution remained nearly constant at ∼16 % between 1971 and 1991 in the Adyl-Su River basin (northern slope of the central Greater Caucasus). Between 1991 and 2006, the supra-glacial debris cover started to increase noticeably, reaching 23 % (∼0.46 % yr −1 ) within 15 years. For the Zopkhito and Laboda glaciers (southern slope of the central Greater Caucasus), supra-glacial debris cover increase was lower in the same period (from 6.2 % to 8.1 %, or ∼0.12 % yr −1 ). We measured the supra-glacial debris cover increase in the same glaciers from 4.9±6.5 % to 6.1±6.4 %, or ∼0.08 % yr −1 , between 1986 and 2000.

5.2 Possible reasons for supra-glacial debris cover changes

We observed a clear increase in supra-glacial debris cover in all investigated regions, which became more pronounced after 2000. Based on our investigation, the upper limit of supra-glacial debris cover migrated up-glacier (Figs. 7, 8; Table S2) as a response to glacier retreat and reduced mass flux, as described by Stokes et al. (2007) and defined as “backwasting” by Benn and Evans (1998). A similar pattern of up-glacier migration has also been described on Tasman Glacier, New Zealand (Kirkbride and Warren, 1999), and on Zmutt Glacier, Swiss Alps (Mölg et al., 2019).

The results presented in this study indicate that the clean-ice area for all selected glaciers decreased from about 93 % to 87 % between 1986 and 2014 (Table 2). This reduction was caused by both glacier retreat and an increase in total supra-glacial debris cover (Table 2; Figs. 3–6, S5). This finding is also consistent with field measurements on Djankuat Glacier (Popovnin et al., 2015).

Glacier thinning and a warming atmosphere can lead to permafrost thawing and slope instability at higher altitudes (Deline et al., 2015). Rock avalanches after 2000 on some glaciers in the Greater Caucasus (particularly in the eastern section) might be one of the reasons why the increase rate was higher during the second time period (2000–2014). For example, a rock–ice avalanche onto Devdoraki Glacier on 17 May 2014 (Tielidze et al., 2019a) caused an area increase in supra-glacial debris cover from 5.9±6.0 % to 19.1±5.6 %, or about 0.95 % yr−1, and a landslide after 2000 onto Suatisi Glacier produced a supra-glacial debris-covered area increase from 2.1±6.1 % to 17.6±5.7 %, or about 1.10 % yr−1 (Fig. S6).

Our investigation shows also that the supra-glacial debris cover increases more quickly on the northern slopes of the Greater Caucasus than on the southern slopes, where higher solar radiation input commonly results in smaller glaciers than on the northern slopes. Furthermore, smaller glaciers on the southern slope exist in high cirques with a much steeper surface. Glaciers on the northern slopes are on average less steep than on the southern slopes mainly because most valley glacier tongues in the north are longer and reach lower altitudes than the south-facing glaciers. This conclusion is supported by Lambrecht et al. (2011), who observed a more rapid increase in supra-glacial debris cover on the northern slopes than the southern slopes.

The highest mean upper limit of the supra-glacial debris cover in the eastern Greater Caucasus can mostly be found at small cirque and simple-valley type glaciers preserved at high altitudes and surrounded by rock walls.

The variation in supra-glacial debris-covered area in the eastern, central and western Greater Caucasus could be mostly caused by climate, lithology and morphological peculiarities of the relief. Some river basins in the eastern Greater Caucasus are built on Jurassic sedimentary rocks, which are characterized by relatively high denudation rates (Gobejishvili et al., 2011; Bochud, 2011) supporting supra-glacial debris cover formation. The relief of the central Greater Caucasus is mainly constructed from Proterozoic and lower Paleozoic plagiogranites, plagiogneisses, quartz diorites, and crystalline slates, which are more resistant and are less prone to the formation of rock avalanches. In addition, the central Greater Caucasus is the highest section of the main range, and glacier tongues are relatively steep and hence less favourable for debris cover accumulation. The western Greater Caucasus is hypsometrically lower, with less steep glaciers. This section is characterized by the highest glacier retreat after the eastern Greater Caucasus. It is therefore possible that glaciers are also rapidly thinning, favouring debris-covered area over the ablation area (Benn et al., 2012; Pratap et al., 2015). The dome of the anticlinorium of the western Greater Caucasus (crest of the main water divide) is built on Proterozoic and Paleozoic plagiogneisses, granites, amphibolites, and crystalline slates. This provides the framework for overall denudation of the high-mountainous relief (over 3000 m; Gobejishvili et al., 2019). Furthermore, this area is characterized by active tectonic and ongoing mountain building (uplifting) processes (Tsereteli et al., 2016), which might be a further reason for increasing supra-glacial debris cover. We note that all these reasons need confirmation by detailed field measurements and could be part of a separate investigation, since there is no accurate geographical pattern which otherwise explains the clear differences of the increase in supra-glacial debris cover.

Our results indicate more than doubling of supra-glacial debris-covered area for Elbrus glaciers from 1986 to 2014, with the highest increase rate between 2000 and 2014 (Figs. 4, 9), although the total estimated uncertainty is comparable to the obtained relative changes. Glaciers on the western slope of Elbrus are affected by avalanches and thus are partially debris-covered (Kutuzov et al., 2019). Glaciers on the eastern slope are characterized by high rates of retreat and great expansion in proglacial lake number and area (Petrakov et al., 2007). The most significant increase in supra-glacial debris cover occurred on the eastern-oriented glaciers of Elbrus, where glaciers are characterized by the highest thinning rates in recent years (Kutuzov et al., 2019). Detailed ground-penetrating radar (GPR) survey helps with more accurately identifying supra-glacial debris cover extent in this area (e.g. GPR measurements by Kutuzov et al., 2019, showed that ∼30 m of ice may be present under the previously considered ice-free area on the eastern slope of Elbrus).

The glaciers in the Greater Caucasus have been retreating continuously since 1960 (Tielidze and Wheate, 2018), suggesting that the shielding effect of increased supra-glacial debris cover at the glacier surface may only partly offset the retreat trend. This is in line with detailed observations of the evolution of Zmutt Glacier, in the Swiss Alps (Mölg et al., 2019), and the fact that mass changes of debris-covered and debris-free glaciers in the Himalaya are similar (e.g. Brun et al., 2019; King et al., 2019). Direct field measurements show that thermal resistance of the <20 cm supra-glacial debris cover for some glaciers (e.g. Djankuat and Zopkhito) in the Greater Caucasus is higher (0.07–0.15 ∘C and 0.05–0.08 ∘C m2 W−1) than in other glacierized regions of the world (e.g. Baltoro, Karakoram: 0.02–0.07 ∘C; Maliy Aktru, Altai: 0.02–0.09 ∘C m2 W−1; Lambrecht et al., 2011), preventing a more rapid retreat. This process is consistent with our observations of the largest debris-covered (Bezingi) and debris-free (Karaugom) glaciers of the Greater Caucasus, where the latter is characterized with higher area shrinkage and terminus retreat. Numerous authors have found similar model results in the Himalaya (e.g. Scherler et al., 2011; Rowan et al., 2015; Jiang et al., 2018).