COMPARISON OF THE VARIABILITY OF SNOW COVER PARAMETERS OF THE HBV MODEL USING LUMPED AND DISTRIBUTED PRECIPITATION INPUTS AND MULTI-BASIN CALIBRATION

Snow cover is a significant source of water supply, mainly in mountainous regions, as snow precipitation fundamentally affects a catchment ́s water balance. The correct simulation of the water balance with rainfall-runoff models is therefore important for the effective management of water resources. Three basic factors may affect the efficiency of hydrological models and the quality of the modelled outputs: The spatial representativeness of the input data, the model ́s structure, and the uncertainties of the model parameters. A comparison of the variability of snow cover parameters and model efficiency of two versions of the HBV model using spatially lumped and distributed precipitation inputs by a multi-basin calibration exercise was performed in this study. Both the lumped and semi-distributed versions of the HBV model were calibrated for discharges, precipitation, and the air temperature on 180 catchments located all over the territory of Austria using data from the period 1991–2000. The analysis focused on the variability of the parameters controlling the snowmelt and the accumulation of the snow components of the two models. The efficiency of the models based on lumped and spatially distributed inputs was compared. The question as to how the catchment ́s mean elevation, and the number of days with an air temperature below zero affects the model ́s performance was targeted, too.


Introduction
Rainfall-runoff models are often applied when solving various water resource problems, e.g., forecasting flood events, estimation of the effects of climate change, simulating extreme discharges, etc. Snow accumulation and snow-melt fundamentally affect the water balance and runoff in river basins. Therefore, the correct modeling of the duration of the snow cover and its properties is an important factor in improving the efficiency of similations of rainfall-runoff models. Three basic factors may affect the efficiency of hydrological models and the quality of the outputs modeled: the spatial representativeness of the input data, the structure of the model, and the uncertainties of the parameters. Model-based approaches are imperfect due to model biases and uncertainties about the input data, too (O'Connell, 1991). In this respect, and based on a review of several studies, Finger et al. (2015), pointed out that analyses of rainfall-runoff model performances in various environments indicated the fact that focusing on a model´s complexity may be less important than the use of proper calibration methods. Kirchner (2006) recommended that in addition to developing better models and better analytic tools, the quality of the data input to models should also get attention. Instead of fitting a pre-defined model and data structure to a catchment via the calibration of parameters, different sources of additional field data could improve the adequacy of representing the dominant hydrological processes in the modelled catchments (McMillan et al., 2011). Respecting the spatial variability of precipitation and snowcover and accounting for other spatially nonhomogeneous basin properties (e.g. soil moisture) could significantly improve the quality of modelled hydrological responses. Khakbaz et al. (2012) discussed issues relating to characterizing the impact of the spatial distribution of rainfall and basin characteristics on runoff generation and the structure of a model. They investigated lumped and distributed calibration strategies and suggested that the performance of a model at an outlet can be improved by using a semi-distributed structure and spatially distributed inputs. A particular problem in the estimation of hydrological model parameters is equifinality, which has been discussed in-depth in a large number of studies (e.g., Freer et al., 1996). As a consequence, many authors have also attempted to constrain such uncertainties in model parameters by using additional data sets for multi-site model calibrations (e.g. Perrin et al., 2001;Hailegeorgis and Alfredsen, 2016;Finger et al., 2015;Knoben et al., 2019).
Research on enhancing the reliability of estimations of snowpack parameters and their contributions to discharges in mountainous regions has also made use of multi-basin data sets. The inclusion of the remote sensing of satellite snow cover images during calibrations has led to the improvement of both snow cover and discharge simulations and the reduction of parameter uncertainties (e.g., Jiang and Wang, 2019;Lopez et al., 2020;Ruelland, 2020). We may account for different climatic conditions in catchments with larger elevation ranges. Various authors have pointed out that not respecting the different climatic conditions in the use of models may lead to uncertainties that could affect the quality of the outputs (e.g., Vaze, 2010;Merz, 2011;Coron, 2012;Saft, 2016;Ceola et al. 2015). Differentiating model inputs by elevation zones may contribute to resolving such problems. One challenge in correctly estimating the properties of snow processes for modelling is the sparse observational network of climatic and hydrological variables in many regions. The preparation of representative model inputs and the resulting model simulations are affected by the uncertainties of the measurements, their spatial representativeness, and the spatial interpolation of the point precipitation measurements. The impact of the spatial representation of the variability of snow cover properties on model simulations therefore continues to receive interest (Finger et al., 2015;Lopez et al., 2020;Ruelland, 2020). Since the snow routine parameters of the conceptual rainfall-runoff models usually cannot be obtained or derived directly from field measurements of the snowpack properties in the climatic stations, they have mainly been estimated by the calibration of the mathematical models. Recently, the streamflowbased model calibrations were extended by remotely sensed snow observations in snow-dominated areas by making use of their increasing spatial resolution and reasonable spatio-temporal coverage. Several studies have shown that incorporating snow observations into the multivariable calibration of a hydrological model could improve streamflow estimates (see Jiang and Wang, 2019). The main objective of this paper is to observe how the lumped and semi-distributed versions of the HBV type TUW conceptual rainfall-runoff model compare when using lumped and spatially distributed climatic inputs in a multi-basin calibration in 180 Austrian catchments. The study closely focuses on a comparison of the parameters of the snow component modeling part of both model versions. The model versions differ mainly in the spatial resolution of the inputs; the semi-distributed version divides each catchment into elevation zones of 200 m, while the lumped model takes every input and output component as a mean value for the whole catchment. The catchments were divided into three groups based on their mean elevation. The model efficiency and snow-related parameters were analyzed separately in catchments with different hypsometric characteristics in flat, hilly, and mountainous catchments. It was attempted to verify if improvement of the model efficiency could be achieved by only using spatially interpolated distributed inputs in a semi-distributed version of a lumped conceptual model (without using distributed parameter calibrations and remotely-sensed snow observations). The analysis also focused on a comparison of the variability of snow cover parameters of two versions of the HBV model using a multi-basin calibration exercise. We hypothesized that the different spatial resolutions of the lumped and semi-distributed models with regard to their input values may lead to observable differences, both in the variability of their performance and parameters, especially in catchments with higher mean elevations in mountainous or alpine regions.

Methods
In this study, the TUW rainfall-runoff model TUW was used in its lumped and semi-distributed versions (Parajka et al., 2007;2009). The model is based on the philosophy of the Swedish HBV model (Bergström, 1995). The lumped version of the TUW model uses the averaged values of the air temperature, precipitation, and potential evapotranspiration as inputs over the whole catchment. The semi-distributed version of the TUW model considers spatially variable inputs over the catchments in 200 m elevation zones. The parameters of the semidistributed model were considered as lumped in this study. Both versions have extensively been used for solving various hydrological problems (see e.g., Sleziak et al., 2016;Parajka et al., 2007;Viglione et al., 2013). The TUW model consists of three sub-models: the snow sub-model, the soil sub-model, and the runoff formation sub-model. Fig. 1 represents the structure of the lumped version of the TUW model. The model has 15 parameters, which are listed in Table 1 together with the recommended ranges of their respective values according to Merz et al., (2011); the same ranges were used in this study. The snow submodel simulates the accumulation of water in a snowpack and inputs water from the melted snow to the catchment. At the centre of interest of this study was the behaviour and variability of the 5 snow routine parameters of both versions of the TUW model. Snow accumulation and snowmelt are controlled by the following parameters:  the snow correction factor (SCF), which represents the uncertainty in the precipitation measurements input in the winter and the large spatial variability of the snow cover; the snowfall is corrected by this corrective snow factor;  the degree-day factor (DDF); a factor influencing the melting of snow;  threshold air temperature (Tr); precipitation above this is considered as rain;  threshold air temperature (Ts); below which, precipitation is considered as snow;  threshold air temperature (Tm); above which, melting in the snowpack takes place.
For calibrating the model in this study, the DEoptim differential evolution algorithm was used (Sleziak et. al., 2017), and the warm-up period was set at one year.
The objective function, which is described in Eq. 1, combined the well-established Nash-Sutcliffe efficiency and the logarithmic Nash-Sutcliffe efficiency (log NSE). The NSE and log NSE coefficients range from -∞ to 1, where 1 indicates a perfect simulation, i.e., an absolute equality between the observed and simulated flows. While the NSE is considered more appropriate for high flows, the log NSE is more appropriate for low flows (Merz et al., 2011).
The following NSE and log NSE formulas were used: where are the simulated mean daily flows, are the observed mean daily flows, is the average of the observed flows.
The objective function (RME) was defined as:

Input data
The calibration of the model was performed on data from 180 catchments, which are distributed over the whole territory of Austria. These data have also been extensively used in previous modeling studies, e.g., by Viglione et al. (2013) and Sleziak et al. (2016). The catchment areas varied from 14.2 km 2 to 6214 km 2 . Before processing the data, quality flags, missing data, etc., were visually inspected. Catchments that were selected which were not affected by an anthropogenic influence, e.g., by dams, canals, or any other artificial runoff regime transformations. The input data (rainfall, runoff, potential evaporation, air temperature) in daily time steps from the period 1.1.1991 to 31.12.2000 were interpolated for the lumped TUW model version from point measurements taken across Austria from 1091 stations by the external drift kriging method (Sleziak et al., 2017). The runoff data were from 180 gauging stations of the Austrian Hydrographical Service. The potential evaporation data were calculated with the Blaney-Criddle method (Parajka et al., 2003). The rainfall and air temperature input data for the semidistributed version of the TUW model were taken from the Spartacus database (Hiebl et al., 2016) and were interpolated into the hypsometric zones by 200 vertical meters. The potential evaporation was calculated with the Blaney-Criddle method in the same hypsometric zones.
In order to separate the effect of the prevailing climatic conditions and the respective runoff regimes in the analysis of the variability of the snow parameters on the results, we clustered the catchments into three groups based on their respective mean elevations: The first group includes catchments where the major Fig. 1. Schematic structure of the lumped version of the TUW model (Sleziak, 2017).
contributor to the runoff is liquid precipitation; this group is labeled the "Lowland" type. The second group includes catchments where a significant part of the runoff is also contributed to by meltwater and is referred to as the "Hilly" type; and the third group represents the "Alpine" type of catchments, where snow and glaciers largely impact the runoff from the catchments. In Fig. 2 we can see the location of the selected catchments, which are clustered into three groups and color-coded as green -Lowland, orange -Hilly, red -Alpine.
The left graph in Fig. 3a shows the median of the number of days with a temperature below zero in the three groups of catchments. The right graph, Fig. 3b, represents the mean elevation of the catchments (726.2 m.a.s.l in the lowland catchments, 1385.7 m.a.s.l in the hilly catchments, and 2212.1 m.a.s.l in the alpine catchments). We can observe that the triangles indicate that a portion of the days with a temperature below zero is directly related to the mean elevation of the catchments; therefore, the separation of the catchments into groups also reflects the differences in the snow regimes.  Illustration of the interdependence between the mean number of days with a temperature below zero and the mean elevation of the catchments divided into the three groups, green -Lowland, orange -Hilly, and red -Alpine catchments.

Results and discussion
Comparing the distribution of the RME efficiency measure of the calibration of both model versions in the groups of catchments (Figs. 4 and 5) shows notable differences in all three groups. The boxplots shown the minimum and maximum value of RME, first and third quartile and mean value represented by black cross. We assessed the model efficiency over three different groups of catchments (i.e., green -Lowland, orange -Hilly, red -Alpine). In the lowland catchments, where the soil-moisture regime is the more dominant runoff generation mechanism, we did not expect a significant contribution of the snow precipitation to the runoff. The RME values were 0.64 in the lumped version of the TUW model and 0.76 in the semi-distributed version (Figs. 4 and 5) in these catchments. A lumped version of the model shows a slightly lower performance in comparison with the semi-distributed version. This means that the semidistributed model also performs better in the lowland catchments, which could be caused by the lower spatial variability of the input data. The Hilly type of catchments with a mean elevation between 1000-2000 m.a.s.l, have seasonal precipitation regime characteristics with the main proportion of liquid precipitation in the summer season and solid precipitation in the winter season. This means that the snow routine of the rainfall runoff model has a stronger influence on the final efficiency of the model´s performance more than in the lowland group of catchments. The median RME values were 0.67 for the lumped version and 0.81 for the semi-distributed version of the TUW model (Figs. 4 and 5). Again, we can observe that the semi-distributed model outperformed the lumped version, probably due to the different spatial resolutions of the inputs into both versions of the model. The third group of catchments with alpine characteristics is the group with snow-dominated catchments, where the melting of accumulated snow precipitation mainly contributes to the catchment´s runoff. In part of the catchments in the alpine group catchments, we may also consider a significant contribution to runoff from glaciers, which represent an important storage of water in high elevation zones. In the alpine group of catchments, the median RME values were 0.51 for the lumped version and 0.88 for the semi-distributed version of the TUW model. Here, we can observe a great difference in the performance between both model versions. The lumped version of the TUW model showed poor performance. The spatial differentiation of the model inputs in the semi-distributed version of the model, which divides catchments into elevation zones, can better reflect the snow regime influenced by the climatic differences between the lowest and highest parts of the catchments. In general, by also comparing Figs. 4 and 5, we can see that the semi-distributed version outperformed the lumped model. We observed that the results of the lumped version of the TUW model showed a poorer performance in catchments with a mean catchment elevation above 2000 m.a.s.l and with alpine climate characteristics (Fig. 4). The semi-distributed version of the TUW model performed better in catchments in the Alpine group compared to those in the Lowland or Hilly groups of catchments (Fig. 5). The main reason for the differences in the quality of the calibration results of both versions of the TUW model may be attributed to the effects of the spatial distribution of the input of the climatic values. Since these values were more pronounced in the Alpine group, it could be expected that the spatially distributed inputs of the winter precipitation improved the model´s performance. In the next step the median values of the snow sub-model parameters were compared in the three groups of catchments (Fig. 6), and the parameter variances in the boxplot charts were compared. The boxplots show the minimum and maximum value of each parameter, first and third quartile and mean value represented by black cross (Fig. 7-11). We can observe differences in the general behavior of the variability of the snow routine parameter among the catchment groups. Fig. 6a shows that the snow  correction factor (SCF) is the highest for both models in the Alpine group; on average, it practically does not strongly correct the winter precipitation amounts in the other two clusters as the boxplot chart (Fig. 7) indicates. This could be expected, and it also shows that the multisite calibration was able to capture this behaviour across all the catchments. The distribution of the values of the degree-day factor (DDF) in Figs. 6b and 8 shows a different pattern: the highest values can be observed for the semidistributed model in the red group, followed by the orange and green catchments. This is consistent with the idea that snowmelt (and snow accumulation, too) has to be most pronounced in an Alpine region, followed by the hilly and lowland catchments. It can be expected in the behaviour of the lumped model that the snow cover may have a greater temporal and spatial variability in the Lowland group. Consequently, the snowmelt can have a shorter duration and be more intensive in some of the lowland catchments, which can explain the distribution and peaking of the median of the DDF parameter there. In the threshold temperature above which precipitation is considered to be liquid (Tr), (Figs. 6c and 9), we can observe consistency in the values for both the semidistributed and lumped versions in the Alpine and Lowland groups, which is to be expected and acceptable. In the Alpine group of catchments the Tr parameter reached almost the same value of +3 o C for all the catchments analysed in this group. The lower median value for the lumped model in the Hilly group of catchments could be connected to the larger variability of the duration and extent of the snowpack at these altitudes, but that would need a detailed catchment-based analysis, which was not performed here. In Figs. 6d and 10, the threshold temperature below which precipitation is snow (Ts) gave similar calibration values for all three groups of catchments and both versions of the TUW model. It is physically acceptable that at lower elevations, the snow accumulation is connected with air temperatures below zero (for majority of chatchments Ts is about -2.5 o C), whereas in high elevations, this temperature can be higher. This could be observed especially in the Alpine group of catchments, where the median value of this parameter reaches a value of about +1 o C. The threshold temperature above which melting starts (Tm) (Figs. 6e and 11), is an important parameter that indicates the start of snow melting and thereby runoff generation. Whereas its similar values in the lumped model´s representation of the spatial variability of the inputs in all three groups are to be expected, the patterns and significant changes in the parameter values in the semi-distributed version of the TUW model are difficult to explain and could be connected to the larger variability of the duration and altitudinal variability of the extent of the snowpack in the orange altitudes (and maybe by the selection of the altitudinal thresholds, too). This would need a detailed catchment-based analysis and maybe a more differentiated subdivision of this group of catchments, which was not performed here.

Fig. 6.
Comparison of the median values of the snow submodel parameters for both model versions. Fig. 7.
Boxplots of the distribution of the SCF parameter in all the catchments divided into the three groups. Green -Lowland, orange -Hilly, and red -Alpine catchments. Fig. 8.
Boxplots of the distribution of the DDF parameter in all the catchments divided into the three groups. Green -Lowland, orange -Hilly, and red -Alpine catchments. Fig. 9.
Boxplots of the distribution of the Tr parameter in all the catchments divided into the three groups. Green -Lowland, orange -Hilly, and red -Alpine catchments. Fig. 10. Boxplots of the distribution of the Ts parameter in all the catchments divided into the three groups. Green -Lowland, orange -Hilly, and red -Alpine catchments. Fig. 11. Boxplots of the distribution of the Tm parameter in all the catchments divided into the three groups. Green -Lowland, orange -Hilly, and red -Alpine catchments.

Conclusion
Snow cover is a significant factor for the supply of water for diverse uses and is an important part of runoff processes, especially in mountainous regions. It is necessary to observe and evaluate how rainfall-runoff models simulate each runoff component in order to ensure reliable simulations that can improve decisions in solving water resources management problems. In this study, we compared the model efficiency of the HBV type TUW rainfall-runoff model in its lumped version and semi-distributed versions. We performed multi-site calibrations of both models for 180 catchments across Austria, which were divided into three groups, according to their respective mean elevation. For the lumped TUW model version, the input data (rainfall, runoff, potential evaporation, air temperature) in daily time steps from the period 1.1.1991 to 31.12.2000 were interpolated from point measurements across Austria (Sleziak et al., 2017) from 1091 stations by the external drift kriging method. The rainfall and air temperature input data for the semidistributed version of the TUW model were taken from the Spartacus database (Hiebl et al., 2016) and were interpolated into the hypsometric zones by 200 vertical meters. We analyzed if and how the mean elevation of the catchments and the spatial variability of the input values can affect both the calibration efficiency and the values of the snow sub-model parameters. The results of the runoff model efficiency showed that the semidistributed version of the model performed better in all the catchments. The efficiency of the lumped version mainly struggled in the group of high altitude catchments. With the variations of the snow sub-model parameters, we can conclude that the overall behavior of the parameter values was physically consistent with the expectations for all the parameters, except for the threshold temperature above which melting starts. Since this can play a huge role in the estimation of the amount of water melted from snow, its behavior in the multi-site calibration requires further analysis, which could improve runoff model efficiency in the semidistributed model version. The spatial differentiation of the model inputs proved to be beneficial and the multisite calibration in the attitudinally grouped catchment clusters led to better insights into the physical consistency and reliability of the snow parameters in the case of the TUW model.