Centenary (1930–2023) climate, and snow cover changes in the Western Alps of Italy. The Ossola valley

In this paper, we study centennial trends of climate and snow cover within the Ossola valley, in the Western Italian Alps. We pursue different tests (Mann Kendall MK, bulk, and sequential/progressive MKprog, Linear Regression, also with change point detection, and moving window average MW) on two datasets, namely (i) dataset1, daily temperature, precipitation, snow depth for 9 stations in the area, during 1930–2018, and (ii) dataset2, snow depth and density, measured twice a month (from February 1st to June 1st) for 47 stations during 2007–2023. We also verify correlation with glacier retreat nearby. In dataset1, we highlight a positive trend for minimum temperature with MK, and Linear Regression. Using MKprog/MW, a negative change of snow cover depth, and duration starting from the late 1980s is found. In dataset2, despite the annual variability in snow cover and 2022–2023 winter drought, we assess the maximum snow water equivalent (SWE) to be delayed with respect to maximum snow depth at high altitude (over a month above 2.700 m a.s.l.), highlighting the effect of settling in decreasing snow depth during spring. We also present a formula linking through Linear Regression the Day of the Year of SWE peak to altitude, relevant to assess the onset of thaw season. Due to the high altitude of the stations, and the paradigmatic nature of the Ossola Valley, hosting Toce River, a main contributor to the Lake Maggiore of Italy, our results are of interest, and can be used as a benchmark for the Italian Alps.

a set of statistical analysis to evaluate trends, discontinuities and variations on series of temperature, precipitation and snow cover. We also assess variation with altitude and date in recent years (2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018)(2019)(2020)(2021)(2022)(2023) of snow depth and density. Being the study area part of the watershed of the Po River, watering a large lowland area for irrigation, drinking water, and hydropower , our results can be used as a benchmark for other studies regarding the Alpine cryosphere, and for water-related assessment in the Po valley.

Study area
Ossola is an alpine valley of the Piedmont region of Italy, covering most of the province of Verbano-Cusio-Ossola (Fig. 1). The province has an area of 2260 km 2 , bordering at the East with Lombardy region and at the North with Swiss Cantons of Valais and Ticino. The valley is shaped by Toce River, springing from Riale Plan, in Formazza Valley (upper part of Ossola Valley), at 1800 m a.s.l., then running for 84 km, collecting water from minor side valleys, before the outlet at Maggiore Lake at 193.5 m a.s.l. The area is of great interest for water resources, given that Ossola Valley nests the Toce River (Ravazzani et al. 2016), a main contributor with Ticino River to the deep subalpine Lake Maggiore, the second largest of Italy after Garda (212 km 2 vs 370 km 2 ). Thence, Ticino, as the outlet of Maggiore Lake, contributes to the Po River.
The region embeds some of the highest and majestic peaks of the Alps, the highest point being the Nordend peak (4609 m a.s.l.) of the Monte Rosa Massif, with an altitude range exceeding 4000 m a.s.l., and high variability in climate conditions. Climate is temperate at Lake Maggiore, with a yearly average air temperature of + 12.4 °C (min + 11.2 °C in 1956, max + 14.3 °C in 2003) during 1956-2006 at the station Verbania-Pallanza (211 m a.s.l.) (Ambrosetti et al. 2006) and extremely cold in Val Formazza, with yearly average temperature of 0.0 °C (+ / − 6.7 °C) at Pian dei Camosci station at 2453 m a.s.l., with a humid temperate and subarctic climate, according to modified Koppen climate classification (Peel et al. 2007). Precipitation is relatively abundant in the Ossola province with 1680 mm on average yearly (estimated with Optimal Interpolation by Loglisci et al. 2012) vs. 1043 mm in the Piedmont region (estimate again with Optimal Interpolation, ARPA 2010).
Referring to the classification by Rau et al. (2005), of the many glaciers within the basin, only two are in the size class 2-5 km 2 . Belvedere and Southern Ohsand have a measured area (in 2014) of 4.51 km 2 , and 2.21 km 2 (Smiraglia and Diolaiuti 2016), and both are strongly retreating. Particularly, Ohsand Glacier is expected to largely shrink in the next 50 years, even with no large temperature increase hereon (Stucchi et al. 2019).

Dataset1
In the following sub-chapters, dataset1 based on daily data is presented together with the applied test statistics and tools for the analysis.

Daily temperature, precipitation and snow depth data (1930-2018)
After the building of some dams for hydropower purposes in Ossola Valley, presently managed by ENEL Green Power company of Italy, weather stations were installed to monitor and forecast water inputs in the artificial lakes. Here, Enel Green Power made available daily records of maximum and minimum temperature T max,min , Precipitation

Temperature data quality assessment
Temperature data were measured with a maximum-minimum mercury thermometer that was reset everyday by the operator. Since no metadata are available, we do not know about the technical characteristics of the instruments, and of changes thereby. As a consequence, we cannot guarantee a priori the accuracy of the dataset (Aguilar et al. 2003), and correct maintenance of the stations, e.g. the absence of erroneous re-positioning that can bias the data (Böhm et al. 2010). We thus preliminarily proceed to verify the homogeneity, and dependability of the data. Internal coherence of the measurements is tested by evaluating the Pearson Correlation Coefficient τ p of annual average T max and T min between the 9 stations. We obtain 9 2 = (9)! (9−2)!(2)! = 36 values of correlation values between all the stations both for T max and T min , whose averages are equal to p,Tmax = 0.27, and p,Tmin 0.43. We also evaluate, as a term of comparison, τ p for annual mean temperature from 48 stations located in Northern Italy of HISTALP dataset (already validated and used, to assess longterm temperature variations in the Alpine area, see Auer et al. 2007). Here, a much higher value is found p,Thistalp = 0.82. Given that the case study stations are located in a much smaller area with respect to HISTALP, and still the data are less correlated, we conclude that our temperature data (or part of them) are somewhat biased. Particularly τ p,av is very low, because maximum temperature measurements, in case of (likely) misplacement of some thermometers, are generally more prone to errors, ca. the double of minimum (Böhm et al. 2010), as they can be affected by the sunlight. Furthermore, like reported by Aguado and Burt (2010), night-time temperature tends to be less variable with respect to daylight, and then minimum temperatures are representative of a larger number of hours. Thus, we neglect maximum temperature measurements in our analysis, and we consider minima only.
We then compare our data against those from the two nearest reference stations displaying continous measurements over the last century. One station is from HISTALP archives, namely Gabiet Lake (Fig. 1, ca. 25 km from nearest station), located in Valle D'Aosta region of Italy, where monthly homogenised mean temperature data are available during 1928-2006. The second one is the Grimsel Hospiz station in Bern Canton of Switzerland ( Fig. 1, ca. 15 km from nearest station) where daily mean temperature is available since 1932 to the present days. We apply the Craddock test (Craddock 1979), which provides a visual assessment of the degree of homogeneity of candidate series, by assessing the cumulative difference between the reference, and the tested series, corrected with a monthly delta to ensure that the mean values of the two series are equal (and the extremes of the series equal to 0). In Fig. 2, we report the output of the test for the two reference stations, highlighting the two series with the smallest discrepancy with respect to the reference series, i.e. Campliccioli and Alpe Cavalli. Because a potential inhomogeneity is detected in the first years of measurements, we neglect measurements before 1949, so largely increasing the correspondence between candidate series and reference ones. Accordingly, test statistics for temperature are applied only to these two stations.

Precipitation correction using fresh snowfall
As well as air temperature, also total precipitation P (mm) is difficult to measure in the high mountains, because of (strong) wind, and low accuracy in the assessment of snow water content by (heated) rain gauges (Kochendorfer et al. 2022), with local studies pointing to an underestimation of snow water content down to − 70% (Mazza and Mercalli 1992). Thus, we consider not reliable the readings of the gauges in the presence of snow and we clean precipitation series by considering null the measurements whenever maximum daily temperature is below 0 °C. In this way, we extract a series of daily rainfall R only, unless for the mixed precipitation events, generally occurring with air temperature near 0°, which here are confused with rain events. Finally, we assess solid precipitation S, i.e. snow contribute, deducing fresh snowfall HN (cm) from the (positive) variations of snow depth (e.g. Bocchiola and Rosso 2007) as follows: Fig. 2 Craddock test of temperature data using as reference stations Lake Gabiet, for the whole period (a) and optimal (b), and Grimsel Hospiz for the whole period (c), and optimal (d). Highlighted in blue and green the two recording stations with best performance where suffix d is for a given day, w is [kg/m 3 ] water density and n is fresh snow density. We use here an average value valid for the Italian Alps, n = 115 kg/m 3 (Valt et al. 2018). After computing liquid and solid precipitation, we can assess total precipitation P (mm) as follows: Once computed precipitation, we also add up as testing variable P days , e.g. the number of days with precipitation larger than 1 mm.

Snow-related variables
In order to examine the seasonality of snow cycle, we consider other additional derived variables, i.e. (i) snow cover duration D0, or number of days in the hydrologic year (September-August), when snow depth is above 0 cm), (ii) snow cover start, and end, SC start and SC end , i.e. the first and last day when snow is present on the ground, (iii) maximum and average snow depth, HS max and HS av [cm], and (iv) average fresh snow, HN [cm]. The variables are also considered at seasonal scale, therefore in Winter (December to February), Spring (March to April), Summer (June to August) and Fall (September to November).

Test statistics
In addition to Long-Term Mean and related standard deviation, to evaluate potential presence of trends for each variable, we use the following tests (see e.g. Bocchiola 2014).
1) Linear Regression, with significance given using p-value, with a level α = 5%. This test allows to highlight (significant) linear trends, which are also pursued using a Change-Point detection procedure, below.
2) Moving Window average MW (De Michele et al. 1998). The moving average is taken over a period of 30 years, as suggested, e.g. by WMO (2017), and then compared with longterm average LT. If MW is not contained within the LT confidence interval with α = 5% the hypothesis of stationarity of the data can be rejected.
3) Mann Kendall MK (Kendall 1975), also in the progressive version (MK prog , used to detect the starting point of a trend). This is a non-parametric test, i.e. it does not take into account any hypothesis about the type (linear/non-linear) of trend, and it is therefore a complement to Linear Regression analysis. 4) Change-Point analysis, including an abrupt change of the mean, MCP, and of slope, LRCP, using Matlab package findchangepts (Killick et al. 2012). This method, frequently used in hydrology (e.g. Wang et al. 2014, Valentin et al. 2018, is based on the hypothesis (1)

Page 8 of 24
Therein, (1) and (2) are the statistical properties of interest (here, mean value, and slope) of a first and a second part of the dataset. Given the probability density function of x * observation P(X * ) , the maximum log-likelihood of the alternate hypothesis H 1 is where * (1) and * (2) are the maximum likelihood estimates of the statistical properties. After locating the change-point k, if the ratio between H 1 and H 0 is above a predetermined threshold, the change of is accepted. Furthermore, to avoid a change-point to be detected in response to (few) anomalous values at the two ends of the series, we decide to neglect the change-points located within 10 years from the start/end of the measurements. Then, an F-Test is pursued, with α = 5%, for the two periods identified by MCP, to verify the presence of a change in the variance of the data, at the change of magnitude point.

Correlation analysis and ice depletion
All the variables analysed here are compared through correlation matrices using the Pearson coefficient τ p , to establish a possible correlation/causality with retreat of glaciers nearby. We use for the purpose the measured annual (outlet length) variation at the front of Muttgletscher (Fig. 1), a small ice body (0.36 km 2 in 2018) located in the Swiss Canton of Valais, at ca. 10 km from the nearest analysed station, the changes of which are registered since 1918 (GLAMOS 2022). The main Italian glaciers in the area, Belvedere and Hosand, are not suitable for this analysis, due both to short/sparse data series (of outlet length), and to their larger size, likely implying a slower response in time of the glaciers (e.g. Bahr et al. 1998;Marta et al. 2021), which makes yearly glacier changes less reactive to local variation in climate (temperature, etc.). Other nearer glaciers (e.g. Corno Glacier, Griessgletscher, Gornergletscher) on the Swiss side were also considered, but due to their larger size and/or lack of snout data therein, they display no meaningful correlation with our dataset.

Dataset2
As done for dataset1, here we illustrate dataset2 and inherent test statistics.

Snow depth and density data 2007-2023
In addition to the daily series, Enel Green Power also provided data of snow depth (43 stations) and density (17 stations), measured manually with snow core samplers (Fig. 1, SM Table 1). The measurements are performed twice a month during February 1 st to June 1 st , i.e. 9 times per year, during 2007-2023. In the 15 stations where we have both snow depth, and density data, we can also evaluate Snow Water Equivalent, i.e. SWE (mm) on the ground as where sw [kg/m 3 ] is the measured snow (pack) density.

Mean and altitudinal gradient
Simple average values of snow depth, snow density and SWE are evaluated for each date during 2007-2023. Due to scarce length of the series, Linear Regression vs time is not recommended here (WMO 2017), but we still pursue a Linear Regression vs altitude to assess change of the considered variables thereby. Particularly, in order to assess end of accumulation period, i.e. start of thawing, we consider through Linear Regression the delay with respect to altitude of the Day of the Year when peak of SWE is reached. Namely, we assess for each station the average through 2007-2023 of the sequential day number of the year (starting day 1 on September 1 st as we consider it as the start of a hydrological year) when assessed SWE is maximum, and then we perform Linear Regression of these day numbers vs. the altitude of the stations. We repeat the same procedure also for HS to verify the correspondence in time of the peak of HS and SWE at different altitude values.

Long-term analysis
We report here in summary Tables 2, 3, 4, 5, 6, and 7 the outcomes of the calculated statistics. For each statistic, unless for Long-Term Mean, we highlight only significant results (a = 5%) as previously reported in Sect. 2.2.5. Due to (bad) quality assessment, as explained, we consider only the values of T min from Alpe Cavalli and Campliccioli, while we neglect the other ones (still reported in italics). Outputs of the MK, and MW tests are given as 1/0, i.e. when non-stationarity is accepted/rejected. In Tables 2, 3, 4, 5, 6, and 7, one can notice how different variables are differently responsive to the considered tests. For instance, HS av and HS max display significant changes according to the Linear Regression, MCP and MK tests, while no significant response is found for the MW, whereas other variables (e.g. S0) change visibly according to the latter test. As a graphical support to Tables 6 and 7, we also report in Fig. 3 the absolute frequency distribution of the years when a change is detected (with MCP, and MK prog ) for main snow variables (i.e. HS av , HS max , D0).  Figure 3 clearly shows that, in 15 cases out of the 32, when changes in one variable are detected, these changes happen in the second half of the 1980s (with one peak in 1987, 9 occurrences). We also show in Fig. 4 the mean annual values of snow depth HS av in our 9 stations. Given the visibly large variability of HS av , a specific year (or a short time window) when abrupt changes happen, they cannot be detected. Local minima can be seen in 1988 that would explain the 3 changing points in 1987 for MCP.

Correlation analysis
In Fig. 5, we show a correlation matrix developed for some of the most significantly correlated variables. The highest correlation is found between snow duration D0, average snow depth HS av , and seasonal cumulated fresh snow HN. Particularly, we find that spring and fall snowfall affects D0 more than winter snowfall, while HS av is affected largely by fall and winter snowfall, possibly remaining on the ground even after spring. Rainfall amount is weakly related against other variables, except for a significant correlation between its value in spring, and fall, i.e. with a delay of one season. Seasonal temperatures are also correlated to each other, but they have apparently weak impact upon snow variables. However, significant correlation of temperature is found against changes of the Muttgletscher front, especially for spring temperature (T spr , − 0.30). Also, slight correlation is seen against the new snow amount in the year HN spr , clearly pointing at an indication of the nourishing/shielding effect of snow for glaciers (e.g. Diolaiuti et al. 2012a, b).
In Table 8, we report mean, Linear Regression coefficient against altitude, and determination coefficient R 2 , for HS, snow on ground density, and the variable SWE assessed. In addition to the mean values, we also report in Fig. 7 the confidence limits as given  Fig. 4 Annual Mean Snow Depth (HS av ) from daily data of 9 stations (in grey). It is also reported the average between the stations (black) and linear trend (dotted black) with coefficient α by standard deviation, and the mean altitude of the measuring stations for each day and variable. HS reaches a mild peak on April 1 st , followed by slow decrease, becoming rapider after May. On the other hand, SWE displays a flatter curve with higher variability during spring, as also noticeable from Fig. 6. The regression coefficients here show increasing values of HS and SWE for higher altitudes. Mean snow density increases from winter to spring, clearly due to snow settling (e.g. Bocchiola and Diolaiuti 2010). However, there is a seemingly little lapse rate against altitude, if not for a quick decrease in May-June. This is likely the result of snow accumulation, larger at higher altitudes and thereby continuing even after spring onset, and snow freeze-thaw cycle, more intense at the lowest altitudes.
Finally, in Fig. 8 we report a (linear) relationship between altitude, the date of maximum snow depth (HS max ) (a), and the date of maximum SWE, i.e. the end date of the accumulation period of snow (b). Particularly, on the y axis, we provide the average Day of the Year (among the sampled ones) when HS/SWE reaches annual peak (HS max , SWE max ), against the altitude of the station in the x axis. We also report a Linear Regression slope β thereby, which displays a delay of peak every 100 m of altitude jump equal to 4.8 for HS and + 7.3 for SWE, with significant determination coefficient R 2 . The two regression lines cross on February 25 th at 1580 m a.s.l., meaning that above that altitude, we assess SWE peak to be delayed of 2.5 days every 100 m with respect to HS peak.

Discussion
In this section, we compare our results with available literature, highlighting potential implications of interest.

Precipitation
Total precipitation, and accordingly number of days with precipitation, does not seem to follow a visible trend according to our Linear Regression test. This result, which could also be due to an approximation in our assessment (i.e. considering constant snow density), is in agreement with several studies for Alpine areas. Among others, Brunetti et al. (2006) analysed long-term precipitation stations in Italy including the ones from the above-mentioned HISTALP dataset, and they found no significant decrease for NW of Italy where dataset1 stations are located. Also, Brugnara and Maugeri (2019) found, in Rosmini Domodossola station in Ossola Valley, only a weak decrease (ca. − 5%) of total precipitation in the last 150 years. Even extending the analysis to the sixteenth century, Casty et al. (2005) were not able to detect a widespread long-term trend, but rather few local anomalies, like we observed here for 5 stations with the MW test. The high interannual variability of precipitation is a main cause of difficult trend detection that also was found to be variable both in space and time in the Alpine area, and again not significant for SW area embedding Ossola Valley (Auer et al. 2007). Even correlation analysis against annual temperature changes, or against global climate indices, such as the North Atlantic Oscillation NAO, gave poor results hitherto (e.g. Schmidli et al. 2002;Bocchiola and Diolaiuti 2010). Moreover, when a nexus between NAO and precipitation was detected, it showed to be unstable, with large changes through time, especially after the 1980s (Colombo et al. 2022).

Snow depth and duration
While yearly mean snow HS av is the result of an interaction between snowfall, compaction and melting throughout the whole year (e.g. Sommerfeld and LaChapelle 1970), HS max is mostly affected by phenomena occurring in the accumulation period, i.e. between (late) October and spring. As such, both variables are frequently used to study cryospheric trends (e.g. Dyrrdal et al. 2013;Marty and Blanchet 2012). As we showed in Tables 2, 3, 4, 5, 6, 7, we find relevant reduction of both variables, either as a trend (Linear Regression, and MK) or as a discontinuity (MKPROG, MCP). In Fig. 3, we also display how most of the changing points here are located in the late 1980s. Indeed, this result seems confirmed by some previous works in the Alpine area. Marty (2008) studied the number of snow days in the Swiss Alps, in 34 stations, and found an abrupt drop after 1988-1989, but with no specific trend, as verified both with MK and Linear Regression. In the Adige valley (NE Italian Alps), Marcolini et al. (2017) observed significantly lower snow depth in 1988-1999 with respect to 1980-1987. Snow cover duration D0 is also quite responsive to climate change (Tables 2, 3, 4, 5, 6, and 7 and Fig. 3) and similarly to HS for several stations changing point is located in the 1980s. Nevertheless, Linear Regression is reliable only for 3 stations below 2000 m a.s.l. and provides a maximum reduction of ca. − 5 days for decade, and still lower than values assessed by Klein et al. (2016) for the Swiss Alps (average − 8.2 days for decade). This inequality worth to be more thoroughly explored, could also be due to the quite different time window considered (from 1930s-1950s to 2018 here vs. 1975-2015 in the Klein study).

Temperature
Linear Trend analysis provides significant results for Campliccioli station, where an increase of + 0.016 °C year −1 during 1949-2018 is detected. This value is slightly lower,  Auer et al. (2007) in SW Alps. We also found a break point for temperature in recent years, in 2008 (mean + 3.4 °C in 1949-2007 vs + 4.4 °C in 2008-2018). Such noticeable jump is coherent with global data, witnessing a temperature increase of + 0.240 °C during the 2000s to 2010s (GISTEMP 2022).

Correlation among variables
Negative, weak correlation was found between temperature and springtime/autumn snowfall, while no correlation was found in winter, when at the high altitude of our stations, respective station altitude. The average HS max and value is represented by the circle size and color and for SWE max it is also reported in the label the effect of warming is still less visible (i.e. temperature remains mostly below 0 °C), as reported, e.g. by Serquet et al. for the Swiss Alps (2011). No significant correlation was found between total precipitation and temperature, with only a small correlation between summer temperature and precipitation, coherently with results from Brunetti et al. (2009) where a similar negative coefficient was found for the SW alpine area, nesting the Piedmont region.
The nexus between glacier retreat and temperature increase is well known and established in literature (Haeberli and Beniston 1998), but the correlation between the two variables at yearly scale cannot always be given for granted. In addition to (winter) precipitation, other phenomena can lead the relationship between glacier retreat and temperature (very) far from linearity. Glaciers' inertia, redistribution of snow by avalanche and wind, debris covering, plus change in ice dynamics, as occurred in the Belvedere Glacier nearby (Haeberli et al. 2002) largely impact annual ice front variation. To partially mitigate these effects, Kendall correlation is considered instead of Pearson, which reduces the impact of eventual outliers due to one of the aforementioned phenomena. Despite limitations in the methodology, Muttgletscher variations are here coherent/correlated with spring and summer temperature, and (less) with seasonal snowfall HN, whereas correlation between (solid) precipitation and ice loss, is well assessed in literature, albeit generally lower than the correlation between the latter and temperature (e.g. Colucci and Guglielmin 2015). Lack of correlation here is most probably due to a larger variability of precipitation with respect to temperature, meaning that values of precipitation assessed in a study area could be less representative of Muttgletscher area, on the Northern side of the Alps. In fact, using the HISTALP dataset we assess correlation between average yearly temperature between NW and SW Alpine region, where Muttgletscher and the study area are located respectively, and we find τ p = 0.93, while on the other side correlation coefficient for annual precipitation for the same regions is much lower τ p = 0.32.

Density and snow water equivalent -dataset2
Maximum snow depth values measured on April 1 st from dataset2 are coherently lower with respect to dataset1 (169 vs 239 cm) since this variable has shown to be quite responsive to climate change (Table 2, picture 3). Because in dataset2, the measurements of snow depth end before total ablation, no conclusions concerning snow cover duration between the two datasets can be set forth. Schöber et al. (2016) analysed the long-term mean (1952-2010) of snow depth, and density for the Tyrolean Alp, at stations between 1400 and 2000 m a.s.l. (vs 1500-2500 m a.s.l. here), and found lower values of snow depth, 0.9 m at April 1 st (vs 1.69 m here), but more similar values of density, with 350 kg m −3 at April 1 st (vs 366 kg m −3 here), and 400 kg m −3 at May 1 st (vs 435 kg m −3 here). While snow density is knowingly conservative (e.g. Bocchiola 2010), higher snow depth values point towards a likely different temperature, and winter precipitation regimes.
Despite the limitation given by the coarse temporal resolution of the dataset and the swinging between wet (e.g. 2009) and dry years (e.g. 2022, 2023), we also notice a lag, increasing with altitude, between maximum snow depth, and maximum SWE, which we assess to be over a month above 2700 m a.s.l. Knowingly in the present literature thaw period is conventionally cast at April 1 st (e.g. Ranzi et al. 1999;Bohr and Aguado 2001), but its dependence upon elevation here depicted was also claimed before (e.g. Bocchiola and Groppelli 2010). The lag of thaw period at high altitude may be likely explained by considerable snowfall contribution even during some springs, maybe camouflaged in terms of snow depth by settling. Another explaining factor can be the refreezing of spring rainfall that then becomes part of the snowpack. This phenomenon was observed, e.g. by Gugerli et al. (2019). They studied a snowfield above Plaine Morte Swiss glacier, at 2690 m a.s.l., and found a similar delay with the peak of SWE in late May, while the measured HS was already on a downward trend since April. A correct assessment of SWE max and of the related date/period can be crucial to predict melt/stream flows in a dry spring/summer, especially in high altitude catchments (Jenicek et al. 2018), whereas better performance in hydrological modeling may be achieved considering spring rainfall as (potential) new SWE when snow cover is present. The relationship here displayed between SWE max date, and altitude seems relevant in this sense.

Conclusions
In this paper, we analyse two high altitude datasets covering Ossola Valley, in the Piedmont region of Italy, namely (i) dataset1, a long-term daily dataset  with measurements of snow depth, precipitation, and temperature, and (ii) dataset2, with more recent (2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018)(2019)(2020)(2021)(2022)(2023) periodic (at fixed dates) data of snow depth and density. The dataset1 dataset, due to the expected evolution in time of measuring tools, and possibly changes in site and setting of the stations, passed through a data quality assessment by correlation analysis and Craddock test. We could highlight changes of behavior of dataset1 data in the second half of the 1980s, when sensibly thinner snow cover than before is detected. Snow depth and density measurements in dataset2 allow us to define an equation linking the date of the annual peak of snow depth and SWE with respect to altitude, with which we also detect a shift between SWE peak and snow depth peak, increasing with altitude. The results we report are of interest as a benchmark in the area, and deserve further investigation, also to assess potential effects thereby upon hydrological cycle in Alpine streams for water resources management under climate change.

Supplementary Information
The online version contains supplementary material available at https:// doi. org/ 10. 1007/ s10584-023-03548-7. are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.