Variations of Lys Glacier (Monte Rosa Massif, Italy) from the Little Ice Age to the Present from Historical and Remote Sensing Datasets

Abstract

Alpine glaciers respond to climate imbalance by adjusting their mass and length. In turn, these changes modify the glacial and periglacial environment, leading to increased supraglacial debris cover, the development of glacial lakes and glacier fragmentation. In this research, we investigated the evolution of Lys Glacier (Monte Rosa Group), by studying length, area and volume changes, and evolution of its supraglacial debris cover and proglacial lakes by means of historical sources and high-resolution aerial and satellite orthophotos. Lys Glacier retreated almost continuously, by nearly 2 km, from its maximum Little Ice Age position. More recently, the glacier lost 11.91% of its area between 1975 and 2014 and underwent fragmentation in 2009. Over the same period, glacier fragmentation and tongue stagnation affected the formation and rapid growth of a series of ice-contact lakes and led to a non-linear debris cover evolution. The glacier was also subjected to strong volume losses, with more than 135 m thinning on the ablation tongue from 1991 to 2014. Analysis of the meteorological records (1927–present) from the closest weather station reveals a considerable increase in average annual temperatures by more than 1°C from the mean of 1971–1989 to the mean of 1990–2017.

1. Introduction

Observations show that the strongest influence of global climate change is recorded in alpine environments and glaciers are the first to be affected by global warming [1]. The response of a glacier to climate imbalance occurs initially through mass loss; eventually, the glacier adjusts to mass changes by changing its geometry, including area and length [2]. These geometric changes are also often accompanied by a shift in the glacier geomorphology, including an increase in the glacier debris cover, which can then decouple the glacier response from the temperature signal, and the appearance of thermokarst features such as kettles and ice-contact lakes [3], with possible impacts on glacier hazards downstream [4]. In fact, debris-covered glaciers and the expansion of supraglacial debris cover on debris-free glaciers are increasingly prominent features of the world’s glaciated catchments in mountain regions including Himalaya [5], Karakoram [6, 7], Andes [8], Alps [9], Southern Alps of New Zealand [10], and Caucasus [11].

Over the past century, Italian glaciers have been shrinking at high rates [12, 13, 14]. Among these, Lys Glacier, located in the Monte Rosa massif (Figure 1), can be considered paradigmatic of the changes affecting alpine glaciers. Its response to increasing temperature and reduced accumulation is in fact similar to that of most glaciers in the Alps and elsewhere [14, 15, 16, 17]: terminus retreat, area reduction and decreasing ice thickness. In addition, Lys Glacier has also recently shown other climate-related changes, including detachment of the main debris-covered tongue from the rest of the glacier body and possible separation of the two branches of the glacier tongue (Figure 1a); variations in supraglacial debris cover (both surface and thickness) on the glacier tongue, thermokarst features (e.g. kettle ponds) and processes (Figure 1b); and calving processes at the terminus where an ice-contact lake developed in the late 1990s (Figure 1c).

Figure 1.

Location of Lys glacier within Italy (red star) and of the weather station used in this study (red dot). The bottom pictures show some typical features found on the glacier tongue in recent years, including (a) separation of two tributaries from the main debris-covered tongue; (b) thermokarst features and processes; (c) calving processes at the ice-contact lake.

In the current context of temperature warming and glacier regression, long-term data reporting changes in glacier properties represent an important asset, as terminus, area and volume fluctuations can provide important information concerning the response of a glacier to climate imbalance. These data are a fundamental input for models that enable reconstructing past climate (from terminus fluctuations [18, 19]) and projecting future glacier changes [20], the availability of meltwater for domestic use and the production of hydroelectric energy [21, 22], as well as the formation of future lakes and potential hazards [23] and the impact of glacier change on tourism [24].

In this study, to describe the evolution of Lys Glacier, a multiple approach was followed. Terminus fluctuations since the early nineteenth century were analyzed using a variety of sources including detailed bulletins with reports of glaciological campaigns; over a more recent period (1975–2014), the area changes of the glacier were estimated by using remote sensing datasets, i.e. satellite and aerial orthophotos, while volume changes were evaluated by comparing a pair of digital elevation models (DEMs) obtained from cartography and satellite images. In addition, these sources permitted us to evaluate the evolution (i.e. surface cover and patterns) of the supraglacial debris cover and of the ice-contact lakes over the same period of observation and gain insights into the geomorphological evolution of the glacier.

2. Study site

Lys Glacier drains the southern flank of the Mont Rosa Group (45°54’ N, 07°50′ E) (Figure 1). The most recent Italian glacier inventory ([14], data from 2005) reports a glacier area of 9.58 km², a south westerly aspect, an elevation range between 2392 and 4323 m a.s.l. and a length of 5.71 km. While Lys Glacier is presently the fourth largest Italian glacier, comparatively few studies have been conducted to investigate its evolution over the past century and recent decades. Strada [25] described terminus fluctuations until the early 1980s; Pelfini et al. [26] estimated the glacier response time using dendrochronology; Rota et al. [27] used DEMs from cartographic sources to calculate volume changes between 1925 and 1994. None of these studies or iconographic sources show evidence of continuous supraglacial debris cover until the late 1980s. Since then, debris cover has been present on the lower sector of the ablation tongue, initially as the continuation of medial moraines at an elevation below 2550 m a.s.l. The medial moraines formed below the glacier icefalls (Figure 2), thanks to the debris supplied by the surrounding deglaciated rock walls (owing to macrogelivation processes, permafrost degradation and structural rock falls). In more recent years, debris has come to cover the entire glacier tongue below 2550 m a.s.l., until its detachment from the upper sector of the glacier, which occurred in 2009 [28]. The occurrence of debris cover is thus a recent phenomenon and probably a consequence of the present deglaciation, which has affected Lys Glacier as well as the other alpine glaciers [1, 15]. Although Lys Glacier cannot be considered a debris-covered glacier even in recent times (as only a small sector of the glacier surface is continuously debris-covered, see [29]), the supraglacial morphologies and the processes affecting the debris-covered glacier tongue are comparable to those of larger debris-covered glaciers in Asia and in the Alps, including thermokarst processes and the formation of cavities and ice-contact lakes [29].

Figure 2.

Comparison of historical photographs of Lys glacier. (a) 1989. The glacier exhibits medial moraines on the distal part of its tongue, which is well developed. (b) 2019. The glacier has retreated above a rock wall in the western sector, leaving a dead ice tongue underneath, with a proglacial lake developing in the depression left by the disappearing tongue. The moraine ridges from the 1980s advance phase are also evident. (photo credits: (a) Willy Monterin-© Archivio Monterin. (b) Fabiano Ventura-© F. Ventura-sulletraccedeighiacciai.com).

3. Materials and methods

3.1 Terminus fluctuations

In this study, the terminus fluctuations of Lys Glacier were retrieved from different sources; early data concerning the nineteenth century are the same used in Strada [25] and references therein, i.e. early cartographic sources (map of the states of Sardinia from 1818 and cadastral maps), photographs and descriptions, particularly those in Monterin [30]; from the early twentieth century, length changes were obtained by extracting the field measurements available from the journals published by the Italian Glaciological Committee (CGI) [31, 32]. Surveys of the glacier terminus positions consist of tape measurements from fixed reference points in the glacier forefield, made at or near the end of the balance year, i.e. in the autumn season, with an estimated accuracy about 0.5 m, which can become worse in case of bad environmental conditions, poorly documented switches to new reference points, measurements taken over very long distances or at points on the snout outside the main flow line, or residual snow patches [12]. While the actual uncertainty from these issues is difficult to quantify, glacier terminus fluctuations remain an important asset to assess global climate and environmental change [19].

We followed recommendations by Citterio et al. [12] to perform quality-checking of the measurements extracted from CGI journals, to minimize the issues related to changing reference points or surveyors over the years. When multiple measurements were available from the same survey, an averaged variation was calculated. In case of missing years in the record, periods up to 5 years where the reference point had remained the same were filled by uniformly distributing the total terminus variation over the missing years.

In case of a gap associated with an undocumented change in a reference point, the gap was not filled, and the first measurement following the gap was considered as the starting point for comparison with the measurement resulting for the following year, as in Citterio et al. [12]. The data record for Lys Glacier from CGI journals is rather uninterrupted and permits analyzing the glacier history in detail over the last century.

3.2 Cartographic analysis

3.2.1 Data sources and preprocessing

Maps, aerial orthophotos and satellite images were analyzed to calculate the area and volume changes of Lys Glacier and changes in the area of the proglacial lakes. The available cartographic sources include large-scale maps from the Val d’Aosta region, produced in 1975 and 1991 at a nominal scale of 1:10000 (Table 1). The maps were available in digital form as rasters, projected in the UTM32N coordinate system based on the ED50 datum; contour lines and elevation points were digitized as shapefiles, and a DEM was produced from the 1991 map with a pixel spacing of 10 m using ArcMap Topo to Raster utility.

Remote sensing data	Date	Resolution (m)	Usage
Regional technical map	1975	10	Glacier outlines
Aerial orthophoto	07/09/1988	0.5	Glacier outlines, debris cover mapping
Regional technical map	1991	10	Glacier outlines, DEM production, volume change
Aerial orthophoto	15/06/1994	0.5	Glacier outlines, debris cover mapping
Aerial orthophoto	31/08/1998	0.5	Glacier and lake outlines, debris cover mapping
Aerial orthophoto	12/09/2003	0.5	Glacier and lake outlines
Aerial orthophoto	04/09/2006	0.5	Glacier and lake outlines, debris cover mapping
Aerial orthophoto	16/07/2012	0.5	Glacier and lake outlines, debris cover mapping
Pleiades PHR1B	01/09/2014	2	Glacier outlines, DEM production, volume change, debris cover mapping

Table 1.

List of datasets used in this study and their usage; the acquisition date (whenever available) is reported as dd/mm/yyyy.

Aerial orthophotos from 1988 to 2012 (Table 1) were obtained from Geoportale Nazionale (http://www.pcn.minambiente.it/mattm/); these have a pixel size of 0.5 m and were downloaded in the UTM32 coordinate system based on the WGS84 datum. Most images were acquired in late summer months (August and September; Table 1), with minimum snow cover, while the orthophotos from 1988 and 2012 were acquired in early summer (June and July) and show evidence of residual snow, which however does not affect the glacier tongue. In addition to the aerial orthophotos, we obtained a stereo pair acquired from the Pleiades satellite constellation from the European Space Agency. The stereo pair was imaged on 1 September 2014 and provided at level 1B processing stage; we further processed it to obtain a DEM and an orthorectified image using NASA’s AMES stereo pipeline (ASP). Eleven ground control points were selected from those available from Val d’Aosta regional authority (http://geonavsct.partout.it/pub/GeoNavITG/monografie.asp) to improve the geolocation accuracy of the DEM through bundle adjustment in the software; the DEM was then produced using the semi-global matching algorithm available in ASP [33], and the raw multispectral images were projected onto the DEM for orthorectification. Both the Pleaides DEM and multispectral image have a final resolution of 2 m.

3.2.2 Calculation of area and volume changes

Based on the available planar data (Table 1), the glacier and proglacial lake outlines were estimated over multiple years: in the 1975 and 1991 technical maps, they were already drawn on the map, while in the other datasets, they were manually digitized as part of this study. The spatial resolution of the orthophotos allowed us to clearly identify the outlines and to distinguish the debris-covered parts of the tongue from the proglacial areas. Before digitization, all available data were reprojected to the UTM32N coordinate system based on the WGS84 datum for consistency. We further evaluated the accuracy of the manual delineation using the buffer method proposed by Paul et al. [34], by allowing the glacier and lake outlines to grow and shrink with a buffer of 2 pixels.

The analysis of elevation and volume changes was based on the comparison of the DEMs from 1991 and 2014. Before performing this comparison, the DEMs were co-registered using the approach developed by Berthier et al. [35] and applied by Fugazza et al. [36]. In this approach, one DEM (‘slave’) is iteratively shifted with respect to a reference DEM (‘master’), to minimize the standard deviation of elevation differences over stable areas located outside the glacier (σ_Δh). We selected the oldest DEM as the reference and resampled both DEMs to a common resolution of 10 m. By applying the co-registration, we obtained a residual σ_Δh of 4.18 m; the average elevation difference over stable areas was 7 m, which was subtracted from the reference DEM. Based on the elevation difference on the glacier surface, we then calculated the volume change as:

$$\mathrm{∆V}={\displaystyle \sum _{\mathrm{k}=1}^{\mathrm{K}}}{\mathrm{d}}_{\mathrm{k}}{\mathrm{A}}_{\mathrm{k}}$$E1

where k = 1 … K are the pixels of the reference surface, A_k is the area of each pixel and d_k is the pixelwise difference between the 1991 and 2014 DEMs.

To express the uncertainty of volume changes, we used the approach described by Fischer et al. [37]: the uncertainty in elevation change (${\mathrm{\sigma }}_{\mathrm{∆h}}$) is not considered as entirely correlated but scaled to account for the effective correlated area A_cor, based on the correlation length (i.e. the distance at which two pixels are correlated). We calculated the correlation length using an empirical semivariogram in R software as 200 m and the effective correlated area as 0.12 km². The uncertainty in volume change ${\mathrm{\sigma }}_{\mathrm{∆h}\_\mathrm{cor}}$ is then expressed as:

$${\mathrm{\sigma }}_{\mathrm{∆h}\_\mathrm{cor}}=\sqrt{{\mathrm{\sigma }}_{\mathrm{∆h}}^2\bullet \frac{{\mathrm{A}}_{\mathrm{cor}}}{5\bullet {\mathrm{A}}_{1991}}}$$E2

where A₁₉₉₁ is the glacier area in 1991. The uncertainty calculation can in our case only be applied to the glacier tongue, as at altitudes above the glacier equilibrium line the presence of snow with little contrast deteriorates DEM reconstruction and leads to larger errors [38], particularly with the oldest DEM from 1991.

3.2.3 Mapping debris cover and its changes

To characterize the evolution of debris cover on Lys Glacier, we adopted the methodology proposed by Azzoni et al. [9] for glaciers in the Ortles-Cevedale Group. Supraglacial debris was mapped by employing a maximum likelihood (ML) supervised classification approach. Four classes were included in the classification: debris, ice, snow and shadow. Although we did not aim at mapping the other surface types, a higher number of classes was chosen to permit an improved mapping of debris cover; the shadow class was also included as some images are affected by topographic shadow; this occurred mostly in the accumulation basins except for the orthophoto from 2003, where shadows occur over the whole glacier; thus, this orthophoto was excluded from the analysis (Table 1). We also excluded the maps from 1975 and 1991, as the information on the presence of debris was not available. For each of the satellite/aerial orthophotos, we independently selected 10 to 15 training areas. To estimate the accuracy of the classification approach, we performed a validation test on the 2006 and 2012 orthophotos, by selecting 100 random points for each test and comparing the results of ML classification against those from a manual classification. The 2 years were chosen because they represent best (2006) and worst (2012) conditions for debris cover mapping, in terms of image quality, presence of snow and shadows. For both years, we calculated overall accuracy (the ratio of the number of correctly classified points to the number of total points), as well as producer’s accuracy (PA) and user’s accuracy (UA), for the different classes. PA denotes the ratio of the number of correctly classified points to the total (manually classified) point for a class and is the complement of omission errors; UA denotes the ratio of the number of correctly classified points to the number of predicted points for a class and is the complement of commission errors.

3.3 Meteorological data

In order to interpret the glaciological data of the study area, we analyzed meteorological variables recorded by the weather station installed in 1927 at Gressoney d’Ejola, 1850 m a.s.l., at a distance of 3.7 km from Lys Glacier (Figure 1). This long dataset allows investigating the recent climate behaviour, which influences the evolution of the glacier. At this weather station, manual observations of air temperature, liquid precipitation, total snow depth and thickness of fresh snow, atmospheric pressure, relative humidity, wind speed and velocity and cloud type and cover were collected three times a day (at 8 am, 2 pm and 7 pm) from 1927 to 2012. The data record is almost uninterrupted except for a short period from 1962 to 1970, when the station was temporarily moved downvalley, at 1730 m a.s.l. in the Orsia village. In 2002, an automatic weather station was installed at the same site, and it replaced the manual weather station after its dismissal in 2012. The automatic weather station measures air temperature, total precipitation (by a heated gauge) and snow height (by an automatic webcam) every hour. The contemporary presence of manual and automatic instruments has allowed ensuring the homogeneity and continuity of the series. Although a weather station also exists at a higher elevation and closer proximity to Lys Glacier (Alpe Courtlys, 1992 m a.s.l., 2.5 km distance), this station was installed in 2001 and does not permit a long-term analysis of climate variables. In this study, we calculated monthly, seasonal and annual averages of air temperatures and monthly, seasonal and annual cumulated liquid/solid precipitation for the analysis of climatological conditions.

4. Results

4.1 Terminus fluctuations

The history of the advances and retreats of Lys Glacier is summarized in Figure 3. We here report the complete record of data but also focus on the most recent period (1975–2017; see box in Figure 3). The general trend is one of retreat, after two short advance pulses during the Little Ice Age (LIA), culminating in 1822, when the glacier reached its peak LIA extension, and 1861. A very negative phase occurred between 1860 and 1882, when the glacier retreated by 941 m. This phase was followed by a short advance pulse at the end of the nineteenth century. The twentieth century saw an initial stability, with a small advance between 1913 and 1922, after which the glacier underwent an almost uninterrupted retreat, with the only exception of the short advance phase 1973–1985. The cumulative retreat since 1812 was 1.55 km, while since 1913 the glacier retreated by 847 m, with an average of −7.99 my^−1.

Figure 3.

Lys glacier terminus fluctuations. The box indicates the observation period analyzed in further detail, showing annual variations as well.

From 1975 to 2017, Lys Glacier retreated for 32 out of the 43 years analyzed, with an average variation of −8.26 my⁻¹. Since the last short-lived advance phase ending in 1985, retreat totalled 443 m (−13.85 my⁻¹). Retreat rates were particularly high between 2003 and 2007, with the most negative variations of −45 m (2007) and − 38 m (2003), while 1991 saw a very short readvance (+8 m), and in 2002 the glacier was considered stable. Retreat rates became lower after 2007, at −7.55 my⁻¹ preluding to the separation of the glacier tongue from the parent glacier.

4.2 Area and volume changes

Over the period of observation, Lys Glacier underwent large changes in area. From 1975 to 1999, little change is evident on the western section of the glacier tongue (Figure 4), which reached the lowest elevations (ca. 2350 m a.s.l.). However, the easternmost part of the glacier tongue had already started retreating by 1999. Far greater changes than in the earlier period occurred from 1999 to 2014, with marked regression both in the western and eastern sections of the glacier tongue (Figure 4a). On the western section, the glacier tongue detached from the upper portion of the glacier at the base of an icefall, whereas on the eastern section the retreat exposed steep rocks underneath the ice. In 2014, the two main branches of the glacier still maintained a small connection at the base of a large rock outcrop (Figure 4a) and in the accumulation basins. Over the years, changes in Lys Glacier area were almost linear in rate, except for the first year, 1975, when the glacier had its largest area (11.3 km²) of the period of observation (Figure 5). Since 1988, the area decreased at a rate of −0.045 km²y⁻¹. In view of the high spatial resolution of the imagery, the uncertainty in the glacier outlines is below 2%, which supports the evidence of the general trend.

Figure 4.

Area changes of Lys glacier and its proglacial lakes over sample years. (a) Lys glacier outlines in 1975, 1999 and 2014. Background is the orthorectified Pleiades satellite image from 2014. (b) the glacier terminus and two proglacial lakes in 2006. (c) the glacier terminus and proglacial lakes in 2014.

Figure 5.

Changes in the area of Lys glacier, its proglacial lake and supraglacial debris cover over the period of observation.

Beside changes in the glacier area, recent years have also seen a transformation of the geomorphological setting of the glacier outwash plain. While the exact date marking the origin of glacier lakes is difficult to establish, by 1999 five small ice-contact lakes had formed on the glacier tongue. Progressive stagnation in the ablation area allowed the lakes to grow and new ones to develop. Thus, from five scattered ponds in 1999 (totalling 1527 m²) a larger overall area is seen together with the coalescence in two larger lakes in 2006 (8730 m² in total). In 2014, three further lakes formed while the existing ones grew in size (Figure 4b,c). The area of glacier lakes increased almost exponentially (summing the area of all lakes) from 1999, reaching 58,560 m² in 2014 (Figure 5). The uncertainty in the lake outlines is larger than that of the glacier outlines, up to 30% in 1999 because of their relatively small size.

The DEM comparison shows a very different picture in the glacier tongue compared to the other areas of the glacier. On the tongue, large negative changes occurred, particularly on the western section, with a maximum of 136.52 m ice loss. On the eastern part of the ablation tongue, the decrease in ice thickness reached instead 41 m. Aside from the glacier tongue, the other areas of the glacier show a more complex pattern, with areas of both positive and negative changes occurring. A loss in glacier thickness is seen particularly in the lower sectors of the western glacier tongue and in a narrow band in the uppermost reaches of the western accumulation basin, as well as on the eastern tongue in its lower easternmost area (Figure 6). A positive difference from 1991 to 2014 can instead be seen particularly on the western accumulation basin, where an apparent increase in thickness is observed between 10 and 20 m, and a maximum of above +75 m was recorded. Considering the entire glacier, the ice loss signal is still predominant, with an average of −4.34 m and a volume change of −47.06 × 10⁶ m³.

Figure 6.

Thickness changes of Lys glacier from 1991 to 2014. Negative values represent loss in ice thickness. The background image is the hillshade produced from the Pleiades DEM acquired in 2014.

4.3 Debris cover evolution

On Lys Glacier, debris initially increased from 1988 to 2000, when it reached almost 1.3 km² on the glacier tongue (Figures 5 and 7a). In 1988, debris was present in narrow medial moraines and more abundantly at the terminus (Figure 2). In 2000, a more homogeneous coverage of the glacier tongue can be seen (Figure 7b), as a mantle of debris appears covering the glacier at the higher elevations at the margins of the western tongue, also suggesting a possible input from the lateral valley walls; conversely, on the eastern tongue and parts of the eastern accumulation basin, coverage is sparse and patchy. Since 2000, the expansion of supraglacial debris appears to have halted as the total area decreased (Figure 5): this effect was probably caused by the shrinking of the stagnating glacier tongue, which was however entirely debris-covered by 2005 (Figure 7c and d). The spread of supraglacial debris appears to have slowed down also at the higher elevations, although limited evidence for increasing coverage is seen for the eastern tongue of the glacier. The relatively high slope and the consequent presence of seracs might limit debris accumulation in those areas.

Figure 7.

Evolution of debris cover for Lys glacier. (a) 1988; (b) 2000; (c) 2005; (d) 2014.

The accuracy of debris cover maps was evaluated separately for 2006 and 2012, as seen in Tables 2 and 3. In both years, PA and UA are very high for debris: both are 100% in 2006, while PA is 89% and UA 97% for 2012. Overall accuracy was 83% in 2006 and 71% in 2012. The limited overall accuracy in both tests is mostly caused by the difficulty in distinguishing between snow and ice (Figure 7b, where no ice was identified), which however is not crucial for the analysis of debris cover. We estimate the accuracy of debris cover mapping in other years to lie between the two values reported for 2006 and 2012; however, it is also possible that the debris cover amount after 2000 was underestimated, owing to the larger presence of snow cover and shadows in areas which were otherwise classified as debris-covered in previous years (compare Figure 7b,c and d).

		Manually classified
		Debris	Ice	Snow	Shadow	Totals
Predicted	Debris	13	0	0	0	13
	Ice	0	21	17	0	38
	Snow	0	0	48	0	48
	Shadow	0	0	0	1	1
	Totals	13	21	65	1	100

Table 2.

Accuracy estimation for debris cover classification based on the aerial orthophoto from 2006.

		Manually classified
		Debris	Ice	Snow	Shadow	Totals
Predicted	Debris	43	0	0	1	44
	Ice	5	8	23	0	36
	Snow	0	0	20	0	20
	Shadow	0	0	0	0	0
	Totals	48	8	43	1	100

Table 3.

Accuracy estimation for debris cover classification based on the aerial orthophoto from 2012.

4.4 Climatological analysis

Due to its high altitude and position, Lys valley experiences cold winters and temperate summers. Heavy rainfall can occur when south humid Mediterranean winds blow and collide against the orographic barrier of the southern slopes of Monte Rosa. Perturbations from the west and northwest are more frequent, but they discharge their rain/snow content mainly on the Mont Blanc and Valais areas, leaving the northeastern extremity of the Aosta Valley almost dry.

From 1928 to 2018 (excluding the period 1962–1970 when the station was temporarily moved downvalley), the mean annual temperature at Gressoney d’Ejola station was +4.4°C ranging from +2.6°C in 1984 to +6.1°C in 2015. Since we cannot focus on the commonly used 30-year reference period (1961–1990), we considered 1952–1961 and 1971–1990. Over this period of observation, the mean annual temperature was +4.0°C, slightly lower than the average of the past 30 years (1989–2018), which was +4.6°C. Climate warming is more evident when summer temperatures are compared: from +11.6°C during 1952–1961/1971–1990 to +12.7°C during 1989–2018, with a mean of +12.2°C over the whole period. The hottest month is July with a mean temperature of +13.2°C followed by August (+12.7°C) and June (+10.7°C). Generally, every 3 years the maximum daily of +25°C is recorded, even if the absolute maximum (+28.2°C) was observed on 4 August 2017 and 11 August 2003; 2003 was the hottest summer on record (average + 15.4°C; Figure 8) with relatively low amounts of liquid precipitation (216.6 mm corresponding to 74% of the mean summer total, 293.6 mm). Conversely, the coldest summer was 1977 (+9.9°C) with very heavy rainfall (428.6 mm corresponding to 146% of the mean summer total).

Figure 8.

Climatological analysis of the Gressoney d’Ejola weather station 1928–2018, including mean summer (JJA—June, July, August) temperature, mean winter (DJF—December, January, February) temperature, cumulative winter precipitation (liquid and solid) and cumulative winter fresh snowfall.

During the winter seasons, the monthly mean temperature is −2.6°C in December, −3.6°C in January and − 2.8°C in February with an absolute minimum up to −25.0°C recorded on 10 February 1986. The mean winter temperature is −3.0°C, and the coldest season was 2009–2010 with an average temperature of −5.6°C (and relatively low precipitation: 174.7 mm of total precipitation, 67.5 cm of snow depth, 181 cm of fresh snow; Figure 8), while the warmest one was 1948–1949 (0.0°C), also characterized by the lowest total precipitation (48.4 mm corresponding to 24% of the total mean winter amount, 201.4 mm; Figure 8), the lowest mean snow depth and cumulative fresh snow (7.7 cm and 52 cm corresponding to 12% of the mean winter amount—63.9 cm—and 28% of the mean winter total fresh snow, 188.7 cm, respectively) and the lowest number of days with snow cover (58 days corresponding to 67% of the mean total winter days, 86.8 days).

Generally, the coldest day is on 5 January (−4.5°C on average). Frost days (T_max < 0°C) generally occur from October to April, even if days with T_min < 0°C can occur even in July. Thaw (T_ave > 0°C) begins at the end of June.

Comparing the two 30-year periods, the mean annual cumulated precipitation (liquid and solid) in 1952–1961/1971–1990 was 1126.9 mm, slightly higher than the amount of 1989–2018 (1090.7 mm). The same results can be observed looking at cumulative fresh snow: from 450.0 cm in 1952–1961/1971–1990 to 385.6 cm in 1989–2018 and from 201.4 cm to 187.5 cm when winter amount are considered. However, the variability remained the same: minima and maxima reached similar values during the two 30-year periods. The maximum total solid precipitation was recorded in winter 1954–1955 (579.7 mm, almost three times the mean value for winter) when very high values of mean snow depth (128.3 cm), cumulative fresh snow (376 cm; Figure 8) and number of days with snow cover (90 days) were observed. In this season, the temperature was equal to the average (−3.0°C).

Considering the whole year (average of 1927–2018), 1134.6 mm of rain and melted snow are distributed in 111 rainy and/or snowy days. Total precipitation is lower than the amount falling in the areas downvalley, which are even more exposed to the southerly humid winds, but much higher than the amount received in the dry Aosta Valley where 500 mm per year are generally recorded owing to its intra-alpine position. From November to April, there is an average 371.7 cm of fresh snowfall each year, almost equally divided between the various months: 61.9 cm per month ranging from 54.9 cm in November to 65.9 in December. Heavy snowfalls occur with Mediterranean temperate humid winds, with the daily maximum of 120 cm recorded on 1 January 1986 and the monthly maximum of 250 cm in April 1989.

5. Discussion

Since the LIA, and particularly over recent times, Lys Glacier has undergone remarkable changes in all the investigated parameters of length, area, volume and debris cover. The terminus fluctuation curve for Lys Glacier is strikingly similar to that of other published curves for Mer de Glace (Mont Blanc region, France) and Unterer Grindelwald (Bernese Oberland, Switzerland), which share some of the longest records of terminus fluctuations for alpine glaciers [39, 40], suggesting that in spite of the distance between these glaciers, the climatic setting with a predominant influence from westerly winds is similar. All three glaciers share two distinct advance phases during the LIA, interrupted by a period of retreat, which according to Vincent et al. [41] was caused by a decrease in winter precipitation; small differences however exist in the timing and magnitude of such advances: the maximum length of Lys Glacier in the past two centuries was reached in 1821, similarly to Mer de Glace (and Rosenlauigletscher, [40]), while the extent of Unterer Grindelwald and most other alpine glaciers peaked around the 1850s [40]. Well documented for all three glaciers is also the rapid retreat following the end of the LIA around 1860, although Mer de Glace also shows a small readvance around 1867, for which no evidence exists for Lys Glacier. Both Lys and Mer de Glace then enter a period of relative stability until the 1930s (although marked by several small advances and retreats); a period of marked retreat follows, attributed to enhanced solar radiation [42] lasting until the short-lived advance phase of the late 1960s/early 1970s, which lasted longer, up to 1995, for Mer de Glace compared to Lys and which is generally observed for mostly alpine glaciers [12], although smaller glaciers tend to have larger readvance periods in the twentieth century, indicating a shorter reaction time [43].

At Gressoney d’Ejola, the late 1970s appear particularly favorable years for glacialism, with cool summers and high amounts of winter snowfall. Since the 1980s, a clear warming trend has emerged for summer temperatures, which are 1.3°C higher from the mean of 1971–1989 to the mean of 1990–2017 (Figure 8) and 1.1°C higher when including 1952–1961/1971–1989 (+12.7°C in 1990–2017 compared to +11.6°C), while no clear signal can be seen in winter temperatures, total and solid precipitation. This is in line with trends previously observed for Italy [44] and high elevation regions [45] and explains the large retreat rates seen since 1985 (−13.85 my⁻¹ in 1985–2017). Considering a longer set of temperature data analyzed for the Alps (1856–1998), the whole twentieth century was characterized by rising temperatures, at a rate of 0.50°C per century (considering summer temperatures [46]), even before the record-breaking decades of the 2000s and 2010s.

As concerns glacier area, the values reported here are always larger than those of the recent Italian glacier inventory (9.58 km² in 2005 [14]), even in 2014, because their outlines did not take into account a small area in the eastern part of the glacier, which was here included for consistency with the 1975 and 1991 outlines drawn on the maps. The rate of change of Lys Glacier found in this study (−0.44 kmy⁻¹ or − 0.4% y⁻¹) is however comparable to that of alpine glaciers, both in Italy and in the other alpine countries [14, 46, 47]. The increase in debris cover is lower than that reported by Azzoni et al. [9] for the Ortles-Cevedale region, where 38 glaciers were reported to have on average a 13.3% higher proportion of their area covered in debris in 2012 (reaching 30%) than 2003, while Lys Glacier went from 7.9% in 1988 to 12.4% in 2000, decreasing again to 8.5% in 2014. These values are also lower than those reported by Shukla et al. [5] for Samudratapu Glacier in Indian Himalayas, where debris cover nearly doubled over less than 3 years and in line with observations of glaciers in Caucasus by Stokes et al. [11], who describe a 3–6% increase in debris cover and a 57% increase in supraglacial and proglacial lake area. The formation of lakes on stagnating debris-covered glacier tongues was also observed by Kirkbride and Warren [48] and is attributed to the presence of ice-cored moraine which prevents meltwater runoff and favors the accumulation of water in depressions left by melting ice. Similarly to Stokes et al. [11], we also found that debris cover has not halted glacier retreat, counter to the evidence that a thick debris cover is known to reduce ablation (see, e.g. [49]). A field campaign conducted in 2006 revealed that debris thickness is generally above 10 cm and up to 60 cm, well above the critical threshold for which the insulating effect prevails on the albedo effect [49]. Thus, while it is possible that mass wasting would have been even higher without debris, it is more likely that glacier retreat occurred owing to the presence of ice-contact lakes and cavities, which enhance melt through backwasting [50].

Unlike retreat rates, mean thickness and volume changes for Lys Glacier are noticeably lower than in similar studies conducted on alpine glaciers: D’Agata et al. [21] report a decrease in ice thickness of −14.91 m for glaciers in Sondrio Province, Central Italian Alps, from 1981 to 2007, corresponding to −0.57 my⁻¹, while we estimated thinning of Lys Glacier to be 0.19 my⁻¹. In their study of all glaciers in the Swiss Alps, Fischer et al. [37] report an area-weighted mass balance of −0.62 m w.e. y⁻¹ from 1980 to 2010, while the geodetic mass balance of Lys Glacier (using a conversion factor of 0.85 accounting for the average density of ice and firn as done by [37]) would be −0.16 m w.e. y⁻¹, which is at the low end of the scale for glaciers analyzed in the study by Fischer et al. [37], although a few Swiss glaciers do share a similar geodetic mass balance. The geodetic mass balance of Lys Glacier would also be lower than that reported by Berthier et al. [51], i.e. –1.05 ± 0.37 m w.e. y⁻¹ for glaciers in the Mont Blanc massif using ASTER, SPOT and Pleiades DEMs between 2000 and 2014. If we exclude the possibility that debris cover has slowed down thinning for reasons stated above, the large differences and the apparent low average elevation changes of Lys Glacier can be explained by (1) the relatively large size of the glacier accumulation basins and the presence of seasonal snow on the Pleiades image (Figure 2a) and (2) interpolation errors, especially in the oldest cartographic DEM, also seen particularly affecting the accumulation basins of the glacier where little change is expected to occur over the years [38], while we observe areas with apparent thickening of 20 m and up to 75 m. Considering the glacier tongue, the glacier thinning is however evident: to compare our findings against those of Rota et al. [27], we selected an area below 2660 m a.s.l., corresponding to 0.7 km². The average ice loss was computed as 62.92 ± 0.81 m from 1991 to 2014, equal to a thinning rate of 2.74 ± 0.03 my⁻¹, and higher than the values reported by Rota et al. [27], i.e. -1.32 my⁻¹ between 1925 and 1953 and − 0.42 my⁻¹ between 1953 and 1994, suggesting an increase in the glacier tongue thinning rates. Our values for the glacier tongue are also in line with the findings of Mölg and Bolch [52] who reported an average elevation change of −67 ± 5.3 my⁻¹ for Zmuttgletscher (Swiss Alps), albeit for a larger area and a longer period, between 1946 and 2005.

6. Conclusion

In this study, we analyzed the evolution of Lys Glacier, one of the largest glaciers in the Italian Alps, by looking at a variety of parameters: terminus fluctuations were studied from historical sources and glaciological bulletins from 1812 to 2017; changes in surface, debris cover and area of supraglacial/proglacial lakes together with volume changes were examined from cartographic and remote sensing datasets from 1975 to 2014. The glacier length variations were found to be similar to those of large glaciers in the Alps such as Mer de Glace and Unterer Grindelwald, indicating a similar climatic setting in spite of the distance of these glaciers; the worst conditions for the glacier development occurred after the end of the Little Ice Age and since 1985 (−443 m from 1985 to 2017) reflecting increasing temperatures as seen from the closest weather station located at Gressoney d’Ejola. Overall, Lys Glacier has retreated by almost 1.6 km since the LIA.

All the other glaciological findings point to a strong glacier reduction, which is interpreted as an evident impact of climate change: the rate of area change was −0.04 km² y⁻¹ since 1988, while glacier volume decreased by −47 × 10⁶ m³ from 1991 to 2014. The glacier debris cover increased from 1988 to 2000, when it covered 12.4% of the glacier area, and then started decreasing again, as a result of glacier shrinking, while the area of the proglacial lakes grew exponentially over the same period. We consider the changes in area and debris cover as highly reliable in view of our accuracy assessment (max 2% error in the glacier outlines and accuracy between 90% and 100% when mapping debris cover), while the uncertainty in volume variations is larger because of the lower quality of the input DEM from 1991.

In view of the present conditions of the glacier, which prevent reaching the glacier tongue, remote sensing remains as the only viable option to investigate the glacier variations in the future, while the detachment of the glacier tongue has further complicated studying terminus fluctuations. To further our understanding of the glacier past conditions, other historical sources should be considered, including pictorial documents to lengthen the record of glacier terminus position and aerial photography from the past century to provide more accurate estimates of volume changes.

Acknowledgments

We acknowledge the DAR—Department of Regional Autonomies—of the Italian presidency of the Council of Ministers for funding this research.

The photographic comparison of Lys Glacier was performed by Fabiano Ventura in the context of the project ‘On the Trail of the Glaciers’— sulletraccedeighiacciai.com.

Pleiades images were obtained from the European Space Agency (project ID 32535). We wish to thank Val d’Aosta region for providing access to the GPS reports and regional technical maps; the meteorological series from Gressoney d’Ejola is available thanks to the efforts and dedication of Umberto and Willy Monterin.

Conflict of interest

The authors declare no conflict of interest.