The HBV model is a bucket-type rainfall-runoff model with a number of routines representing the main components of the terrestrial part of the water cycle, i.e. snow routine, soil routine, groundwater (response) routine, and routing function. In this study, we focused solely on the snow routine of the model. We used the HBV-light software, which follows the general structure of other implementations of the HBV model and includes some additional functionalities such as Monte Carlo runs and a genetic algorithm for automated optimisation (Seibert and Vis, 2012). Henceforth we use the term “HBV model” when referring to our simulations using the default HBV-light software, that is, the HBV model with the snow routine as described in Lindström et al. (1997) and Seibert and Vis (2012).

2.2 Proposed modifications to individual components of the snow routine

Here we review the different components of the snow routine structure of the HBV model as well as functions that are directly related to this routine (e.g. input data correction with elevation) and describe the proposed modifications to each component. Each of these alternative structures requires one to three additional parameters (Table 1).

Table 1Description of the proposed modifications to the snow routine of the HBV model. The default component structures of the HBV model are marked with a ¹ symbol. The components marked with a ² are not part of the snow routine but were included in the analysis due to their significant impact on it.

Download Print Version | Download XLSX

2.2.1 Temperature and precipitation lapse rates

When different elevation zones are used, the temperature for each zone is generally computed from some catchment-average value and a lapse rate parameter. In HBV, a constant temperature lapse rate is usually used. Alternatively, if the available data allow, it is also possible to provide an estimation of the daily temperature lapse rate. However, if no data on the altitude dependence of temperature are available, setting a constant value throughout the year might be an oversimplification. Indeed, in an experimental study on several locations across the Alps, Rolland (2002) found that the seasonal variability of the temperature lapse rate follows approximately a sine curve with a minimum around the winter solstice. Following these findings, we implemented a seasonally variable computation of the temperature lapse rate using a sine function (Eq. 5). This way, the temperature lapse rate for a given day of the year, Γ_n (^∘C / 100 m) (where n is the day of the year, a sequential day number starting with 1 on 1 January), depends on two parameters, namely the annual temperature lapse rate average, Γ₀ (^∘C / 100 m), and amplitude, Γ_i (^∘C / 100 m).

$$\begin{array}{}\text{(5)}& {\mathrm{\Gamma }}_n={\mathrm{\Gamma }}_{\mathrm{0}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\mathrm{\Gamma }}_i\sin {\displaystyle \frac{\mathrm{2}\mathit{\pi }(n-\mathrm{81})}{\mathrm{365}}}\end{array}$$

Precipitation lapse rates cannot be related to seasonal or other types of systematic variations as they are strongly dependent on the synoptic meteorological conditions and therefore highly variable. Consequently, we decided to keep the default approach in the HBV model which consists in calibrating the model using a constant precipitation lapse rate parameter.

2.2.2 Precipitation phase partition

The determination of the precipitation phase is a crucial step as it controls whether water accumulates in the snowpack or contributes directly to recharge and runoff. In the HBV model, the distinction between rainfall and snowfall is based on the assumption that precipitation falls either as rain or as snow, depending on a threshold temperature parameter. However, in reality, this transition is less sharp, as both rain and snow can coincide (Dai, 2008; Magnusson et al., 2014; Sims and Liu, 2015). Additionally, depending on other factors such as humidity and atmosphere stratification, the shift from rain to snow can occur at different temperatures. Therefore, the single threshold temperature may not adequately represent the snow accumulation, especially in areas or periods with temperatures close to 0 ^∘C. Different formulations have been proposed to describe the snow fraction of precipitation, S (–), as a function of temperature (Froidurot et al., 2014; Magnusson et al., 2014; Viviroli et al., 2007). In this study, we considered three different formulations to calculate the snowfall fraction of precipitation (Eqs. 6 and 7 respectively): (i) a linear function (Eq. 8), (ii) a sine function (Eq. 9), and (iii) an exponential function (Eq. 10). Both the T_A (^∘C) and M_P (^∘C) parameters control the range of temperatures for mixed precipitation.

$$\begin{array}{}\text{(6)}& {\displaystyle }& {\displaystyle }{P}_{\mathrm{S}}=P\cdot S\cdot C_{\mathrm{SF}},\text{(7)}& {\displaystyle }& {\displaystyle }{P}_{\mathrm{R}}=P\cdot (\mathrm{1}-S),\text{(8)}& {\displaystyle }& {\displaystyle }S=\left\{\begin{array}{ll}\mathrm{1},& \text{ if }T\le T_{\mathrm{T}}-\frac{T_{\mathrm{A}}}{\mathrm{2}}\\ \frac{\mathrm{1}}{\mathrm{2}}+\frac{T_{\mathrm{T}}-T}{T_{\mathrm{A}}},& \text{ if }{T}_{\mathrm{T}}-\frac{T_{\mathrm{A}}}{\mathrm{2}}<T\le T_{\mathrm{T}}+\frac{T_{\mathrm{A}}}{\mathrm{2}}\\ \mathrm{0},& \text{ if }T>T_{\mathrm{T}}+\frac{T_{\mathrm{A}}}{\mathrm{2}}\end{array}\right.,\text{(9)}& {\displaystyle }& {\displaystyle }S=\left\{\begin{array}{ll}\mathrm{1},& \text{ if }T\le T_{\mathrm{T}}-\frac{T_{\mathrm{A}}}{\mathrm{2}}\\ \frac{\mathrm{1}}{\mathrm{2}}-\frac{\mathrm{1}}{\mathrm{2}}\sin \left(\mathit{\pi }\frac{T_{\mathrm{T}}-T}{T_{\mathrm{A}}}\right),& \text{ if }{T}_{\mathrm{T}}-\frac{T_{\mathrm{A}}}{\mathrm{2}}<T\le T_{\mathrm{T}}+\frac{T_{\mathrm{A}}}{\mathrm{2}}\\ \mathrm{0},& \text{ if }T>T_{\mathrm{T}}+\frac{T_{\mathrm{A}}}{\mathrm{2}}\end{array}\right.,\text{(10)}& {\displaystyle }& {\displaystyle }S={\displaystyle \frac{\mathrm{1}}{\mathrm{1}+e^{\frac{T-T_{\mathrm{T}}}{M_{\mathrm{P}}}}}}.\end{array}$$

2.2.3 Snowmelt threshold temperature

In addition to determining the precipitation phase, a temperature threshold parameter is also needed to determine the onset of snowmelt. The most straightforward approach, used in the HBV model, is to use the same threshold temperature parameter for both snowfall and snowmelt. However, as these two transitions are related to different processes happening at different environmental conditions, a single parameter might not adequately describe both transitions. A more realistic approach would be to consider two separate parameters for these processes: a threshold temperature parameter for precipitation phase partitioning, T_P (^∘C), and another one for snowmelt and refreezing processes, T_M (^∘C) (Debele et al., 2010).

2.2.4 Degree-day factor

The degree-day factor is an empirical factor that relates the rate of snowmelt to air temperature (Ohmura, 2001). In the HBV model, a simple proportionality coefficient to estimate the magnitude of the snowmelt is used. This coefficient, multiplied by a constant (usually set to 0.05 in HBV), is also used to compute refreezing rates. Nevertheless, while the degree-day factor is often assumed to be constant over time, seasonal changes in snow albedo and solar inclination point to temporal variations of the degree-day factor as well. While some models use monthly values for this parameter (Quick and Pipes, 1977), a more elegant but still simple way to represent this variability is to consider a seasonally variable degree-day factor following a sine function defined by a yearly average degree-day factor parameter, C₀ (mm Δt^−1 ∘C⁻¹), and an amplitude parameter, C_0,a (mm Δt^−1 ∘C⁻¹), defining the amplitude of the seasonal variation (Eq. 11) (Braun and Renner, 1992; Hottelet et al., 1994). By establishing a seasonally variable degree-day factor instead of a constant value for this parameter, potential snowmelt rates become smaller during the winter months and increase during spring and summer (if there is any snow left).

$$\begin{array}{}\text{(11)}& C_{\mathrm{0},n}=C_{\mathrm{0}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{C}_{\mathrm{0},\mathrm{a}}\sin {\displaystyle \frac{\mathrm{2}\mathit{\pi }(n-\mathrm{81})}{\mathrm{365}}}\end{array}$$

2.2.5 Snowmelt and refreezing

All liquid water produced by snowmelt does not leave the snowpack directly, as a certain amount of liquid water can be stored in the snow, thus delaying the outflow of water from the snowpack. The liquid water stored in the snowpack can also refreeze if temperatures decrease below the freezing point. In the HBV model, both the storage of liquid water and refreezing processes are considered. However, since the magnitude of refreezing meltwater is generally tiny compared to other fluxes, some models disregard this process entirely to reduce model complexity (Magnusson et al., 2014). Here we follow the approach by Magnusson et al. (2014) which, besides disregarding the refreezing process, describes the snowmelt magnitude using an exponential function (Eq. 12). This formulation of snowmelt is somewhat more detailed than the one used in HBV and requires the use of an additional parameter to control for the smoothness of the snowmelt transition, M_M (^∘C). In contrast to the formulation used in the standard HBV model, snowmelt occurs even below the freezing point, but at negligible amounts. The impact of increasing temperature on snowmelt is higher for this formulation compared to HBV.

$$\begin{array}{}\text{(12)}& M=C_{\mathrm{0}}\cdot M_{\mathrm{M}}\left[{\displaystyle \frac{T-T_{\mathrm{T}}}{M_{\mathrm{M}}}}+\ln \left(\mathrm{1}+e^{-\frac{T-T_{\mathrm{T}}}{M_{\mathrm{M}}}}\right)\right]\end{array}$$

2.3 Study domain and data

We selected two sets of mountainous catchments located at different countries within central Europe to test the proposed modifications to the individual components of the snow routine of the HBV model (Table 2, Fig. 1). The first set, composed of Swiss catchments, contains catchments ranging from high-altitude, steep catchments in the central Alps to low-altitude catchments in the Pre-Alps and Jura mountains with gentler topography. The second set, composed of Czech catchments, is representative of mountain catchments at lower elevations compared to the Swiss catchments.

Figure 1Geographical locations of the catchments used in this study. We used a total of 54 catchments, 22 located in Switzerland and 32 in Czechia.

Table 2Relevant physical characteristics of the catchments included in the study. Each catchment is given an identification code in the following way: country (CH – Switzerland, CZ – Czechia), geographical location (Switzerland: 100 – Jura and Swiss Plateau, 200 – Central Alps, 300 – Southern Alps; Czechia: 100 – Bohemian Forest, 200 – Western Sudetes, 300 – Central Sudetes, 400 – Carpathians) and a sequential number for increasingly snow-dominated catchments within each geographical setting. The official hydrometric station IDs from FOEN and CHMI are also provided (four- and six-digit numbers in column station ID).

Download Print Version | Download XLSX

2.3.1 Switzerland

We selected 22 catchments in Switzerland covering a wide range of elevations and areas in the three main hydro-geographical domains of the country, i.e. the Jura and Swiss Plateau, the Central Alps, and the Southern Alps (Weingartner and Aschwanden, 1989). No catchments with significant karst or glacierised areas, as well as catchments with substantial human influence on runoff, were selected for this study. This decision allowed us to observe the signal of snow processes, without including noise or added complexity from other processes, but limited the number of catchments in high altitudes, which are the ones with the largest snowmelt contribution to runoff. The resulting set of catchments had mean elevations between 600 and 2800 m a.s.l. with elevation gradients of up to 2000 m and catchment areas between 3 and 500 km² (Fig. 2). There was a considerable variability in the yearly snowmelt contribution to runoff, ranging from 5 % to 60 % as the catchments ranged from pluvial to glacio-nival regimes.

Figure 2Distribution of the catchments used in this study in terms of area (x axis), mean elevation and range (y axis), and relative snowmelt contribution to runoff (ranging form 4% to 62%; indicated by circle size). The catchments are coloured according to their respective geographical domain: blue (Switzerland) and orange (Czechia).

Download

We obtained the necessary meteorological data for running the HBV model from the Swiss Federal Office of Meteorology and Climatology (MeteoSwiss). More specifically, we used pre-processed gridded data products to obtain catchment-average precipitation (Frei et al., 2006; Frei and Schär, 1998) and temperature (Frei, 2014). These gridded data products are available from 1961 and have a daily temporal resolution and a spatial resolution of 1.25^∘ minutes covering the entire country.

We used both stream runoff and snow water equivalent data for model calibration and validation. We obtained daily stream runoff data from the Swiss Federal Office for the Environment (FOEN). Regarding snow water equivalent, we used 18 years of gridded daily snow water equivalent data at 1 km² resolution derived from a temperature-index snow model with integrated three-dimensional sequential assimilation of observed snow data from 338 stations of the snow monitoring networks of MeteoSwiss and the Swiss Institute for Snow and Avalanche Research (SLF) (Griessinger et al., 2016; Magnusson et al., 2014). Even if using a temperature-index model for both the HBV model and the estimation of the snow water equivalent validation data may introduce some bias to the results, the data assimilation and error correction methods used in the estimation of snow water equivalent make this methodology especially robust (Magnusson et al., 2014). Finally, we obtained the catchment areas and topography from a digital elevation model with a resolution of 25 m from the Swiss Federal Office of Topography (swisstopo).

2.3.2 Czechia

The second set of catchments was composed of Czech catchments and includes 32 mountain catchments with catchment areas ranging from 3 to 383 km² (Fig. 2). As for Switzerland, we selected near-natural catchments with no major human influences such as big dams or water transfers. The selected catchments were located at lower elevations than most of the selected Swiss catchments. Additionally, they were located in the transient zone between oceanic and continental climate, with lower mean annual precipitation than the Swiss catchments. The mean annual snow water equivalent peak for the period 1980–2014 ranged from 35 to 742 mm depending on catchment elevation, resulting in 13 % to 39 % of the annual runoff coming from spring snowmelt.

We obtained daily precipitation, daily mean air temperature, and daily mean runoff time series from the Czech Hydrometeorological Institute (CHMI). Additionally, we obtained weekly snow water equivalent data from CHMI (measured each Monday at 07:00 CET). Since no gridded precipitation, air temperature, or snow water equivalent data are available for Czechia, station data were used for HBV model parametrization. We used stations located within the individual catchments when available. If no such station was available, we selected the nearest station representing similar conditions to the target catchment (e.g. stations situated at a similar elevation). Finally, we used a digital elevation model with a vertical resolution of 5 m from the Czech Office for Surveying, Mapping and Cadastre to obtain catchment areas and elevation distributions.

The large number of catchments and model variations considered in this study made it challenging to grasp any details when looking at the entire dataset. For this reason, we first present the results for a single catchment to explore the implications of individual model modifications and illustrate the general trends observed across the study domain. For this purpose we selected the Allenbach catchment at Adelboden (CH-204), one of the high-altitude, snow-dominated catchments in the set, as the sample catchment. Thereafter, we progressively add more elements to the analysis of the results. Additionally, even if we calibrated (and validated) the model for both periods defined in the split-sample test, here we only present the results for the calibration effort in period 1 and corresponding model validation in period 2, as they are representative of the entire analysis. A comprehensive list including calibration and validation model performance values for both objective functions and all catchments included in this study can be found in Appendix A.

The calibration performance of the standard HBV model for the Allenbach catchment was satisfactory for both objective functions (model performance values of ∼0.90) but, still, some modifications led to increased model performances. Amongst the different changes in single components of the snow routine structure of the HBV model that we evaluated in this study, using a seasonally varying degree-day factor (C_0,s) had the most substantial impact on the model performance to represent snow water equivalent followed by, to a lesser extent, stream runoff (Fig. 3). Apart from this modification, only the use of an exponential function to define the precipitation partition between rain and snow (ΔP_e) produced significant changes in the model performance against both objective functions. In this case, however, this modification impacted the model performance in opposite ways, leading to decreased model performance for the calibration against stream runoff. Model uncertainty, as given by the performance ranges obtained when aggregating the different calibration efforts, was small when compared to the performance differences between the different model structures.

Figure 3Model calibration performance for the 10 calibration efforts against the two objective functions (a: snow water equivalent; b: logarithmic stream runoff) for each of the modifications to individual components of the snow routine of the HBV model for the Allenbach catchment at Adelboden (CH-204). The modifications include a seasonally variable temperature lapse rate (Γ_s), a linear, sinusoidal, and exponential function for the precipitation phase partition (ΔP_l, ΔP_s, and ΔP_e respectively), different thresholds for precipitation phase and snowmelt (T_P,M), a seasonally variable degree-day factor (C_0,s), and an exponential snowmelt with no refreezing (M_e). The median performance of HBV is represented with an orange horizontal line.

Download

Looking at a sample year within the calibration period, we can get a grasp on how comparable the simulated values of snow water equivalent and stream runoff (including the model uncertainty) are to the observed values, both for model calibration and model validation (Fig. 4). While capturing the general evolution of the snowpack, the HBV model tended to underestimate the snow water equivalent amounts, except for the spring snowmelt period. The model alternative using a seasonal degree-day factor (C_0,s), which had proven to be the best possible model structure modification for model calibration against snow water equivalent for this catchment, exhibited the same overall behaviour but was more accurate and precise than the HBV model. Regarding the calibration against stream runoff, both model alternatives performed well for low-flow periods, but they missed or underestimated some of the peaks. Model uncertainty was comparable for both model alternatives and was not significant when compared to the simulated values. Model results from the same sample year for all the catchments included in this study can be found in the Supplement.

Figure 4Example time series (October 2003–September 2004) from the Allenbach catchment at Adelboden (CH-204). (a) Daily mean air temperature and total precipitation. (b, c) Model calibration results for period 1. (d, e) Model validation results (based on model calibration on period 2). The model calibration and validation are further subdivided into (a) catchment-average observed (grey line) and simulated snow water equivalent (HBV in blue and the model structure modification including a seasonally varying degree-day factor, C_0,s, in orange) and (d, e) observed (grey line) and simulated stream runoff (HBV in blue and the model structure modification including a seasonal degree-day factor in orange). The grey field represents the period used when calibrating the model against the logarithmic stream runoff. The uncertainty fields for model simulation cover the 10th–90th percentile range, while the solid line represents the median value.

Download

Looking at the entire set of catchments, the impact of the different model structure modifications on model calibration performance was generally more pronounced for R_W than for R_ln (Q) across all catchments (Fig. 5). For most catchments, the largest model performance improvements when calibrating against snow water equivalent were achieved by using a seasonally variable degree-day factor (C_0,s). Using different thresholds for precipitation phase partition and snowmelt (T_P,M) and using an exponential function for precipitation phase partition (ΔP_e) also conveyed a significant improvement for some of the catchments. Nevertheless, the latter modification performed almost equally to the HBV model when calibrating against stream runoff, and even slightly worse for some catchments. Using an exponential function to define the precipitation partition between rain and snow consistently penalised the model performance when calibrating the HBV model against stream runoff, whereas using an exponential function for snowmelt (M_e) was the best alternative when calibrating the model against this objective function. Overall, most modifications conveyed slight model performance improvements concerning snow water equivalent simulations for most of the catchments in the dataset. Nevertheless, the modifications to the precipitation phase partition tended to penalise most Czech catchments when calibrating against snow water equivalent. We did not observe any significant connection between model performance and catchment characteristics such as mean catchment elevation, catchment area, or yearly snowmelt contribution to runoff.

Figure 5Median relative model calibration performance for alternative HBV model structures, including modifications to single components of the snow routine with respect to HBV. The modifications include a seasonally variable temperature lapse rate (Γ_s), a linear, sinusoidal, and exponential function for the precipitation phase partition (ΔP_l, ΔP_s, and ΔP_e respectively), different thresholds for precipitation phase and snowmelt (T_P,M), a seasonally variable degree-day factor (C_0,s), and an exponential snowmelt with no refreezing (M_e). (a) Model calibration against snow water equivalent; (b) model calibration against logarithmic stream runoff. The catchments are ordered by mean yearly snowmelt contribution to runoff in downwards increasing order.

Download

While some modifications had a clear and consistent impact on model calibration performance in all catchments, most of them presented a less pronounced either positive or negative impact, depending on the catchment, making it difficult to evaluate which model structures were more suitable than others (including the default HBV structure) for most of the catchments. Additionally, to better understand the usefulness of the different modifications in real applications, we need to take into account which of the model structures performed best for the validation period as well (Fig. 6). As already observed in Fig. 5, using a seasonal degree-day factor (C_0,s) was the best modification for calibrating the model against snow water equivalent for the vast majority of the catchments. Nevertheless, this modification ranked relatively low when validating the model against the same objective function. Looking at stream runoff, using an exponential function for snowmelt simulation while disregarding the refreezing process (M_e) was the best-ranking modification for both model calibration and validation, while the HBV model ranked higher than several of the considered modifications. Using an exponential function to define the precipitation partition between rain and snow (ΔP_e) was the worst alternative for calibrating the model against stream runoff. The diagonal pattern from the top left to the bottom right observed for model calibration indicates that modifications tended to have the same rank for most catchments (note that the ranking of modifications is different when looking at snow water equivalent with respect to stream runoff). Such a pattern was not present for model validation, suggesting that, in that case, there were no model structures significantly more suitable than others.

Figure 6Rank matrices for each of the model simulation scenarios. (a, c) Model calibration (a) and validation (c) against snow water equivalent; (b, d) model calibration (b) and validation (d) against logarithmic stream runoff. Each rank matrix shows the rank distribution of each modification to single components of the snow structure of the HBV model for all the catchments included in this study so that each column adds to 100 %. The modifications are ordered from highest to lowest average ranking (left to right) and include a seasonally variable temperature lapse rate (Γ_s), a linear, sinusoidal, and exponential function for the precipitation phase partition (ΔP_l, ΔP_s, and ΔP_e respectively), different thresholds for precipitation phase and snowmelt (T_P,M), a seasonally variable degree-day factor (C_0,s), and an exponential snowmelt with no refreezing (M_e). The standard HBV model structure is highlighted with a white vertical line.

Download

Even if some alternative model structures clearly improved the calibration performance of the model, most structures had a limited impact on model performance. This is in part because, to this point, we only tested model structures containing a single modification with respect to the HBV model. We next explored whether the same trends persisted when including further elements to the model by using an increasing amount of model structure modifications simultaneously. In total we tested 64 different model structures (Table 3). Using the maximum possible number of simultaneous modifications (Eq. 5) to the model structure would result in the use of up to nine additional parameters. This would lead to a clear overparameterisation of the model (the default snow routine structure of HBV contains five parameters), but we included this alternative to provide a complete analysis of all the available alternatives.

Table 3Number of model structure alternatives containing a given number of snow routine modifications.

Download Print Version | Download XLSX

Figure 7 shows the median model performance for each of the 64 possible model structure alternatives for all catchments relative to the standard HBV model performance, sorted by the number of components being modified. When calibrating the model against snow water equivalent, model performance clearly increased for all of the model structure alternatives. The impact was more modest for model validation, with a significant percentage of alternative structures performing worse than HBV. Regarding model calibration against stream runoff, the effect of an increasing number of components being modified was limited but mostly positive. The range of model performance values was also significantly smaller than when looking at snow water equivalent. This relates to the fact that, for most catchments, the snow routine has a limited weight over the entire HBV model. For model validation we observed a similar trend but with broader model performance ranges. Also, the fact that performance variability varied significantly with the number of components being modified was in part due to the differences in the number of model structure alternatives for each of them, being larger for the number of modifications which included the largest number of model structure alternatives.

Figure 7Median relative model performance with respect to HBV across all catchments for the 64 considered snow routine structures sorted by an increasing number of modifications. Model simulations against snow water equivalent are presented in (a, c), while those against logarithmic stream runoff are shown in (b, d). Model calibration is presented in (a, b) and model validation in (c, d). The dashed line represents the median value across all catchments, while the grey fields represent the minimum to maximum (light grey) and 25th to 75th percentiles (dark grey). The relative model performance of HBV is highlighted with a solid orange line.

Download

In general, we observed an increase in model performance for all cases. However, except for model calibration against snow water equivalent, there was no clear indication that a model with a more detailed formulation of the snow processes would lead to significantly improved model performance. Indeed, the range of performances among the different model structures was larger than the median net increase. This might be an indication that choosing the right modifications (or combinations of modifications) is more relevant than significantly increasing model detail. This way, we attempted to determine whether some specific modifications conveyed a model performance gain across all model structures in which they were included for both model calibration and validation against the two objective functions. To this purpose, we ranked all model structures for model calibration and validation against both objective functions and visualised the cumulative distribution of each of the individual model modifications (Fig. 8). Some of the patterns observed here resemble those that we already observed for single modifications only (Fig. 6). For instance, all top-ranking model structures included a seasonally variable degree-day factor (C_0,s) for model calibration against snow water equivalent. Similarly, all bottom-ranking model structures used an exponential function for precipitation phase partition (ΔP_e) for model calibration against stream runoff. Besides these familiar patterns, other patterns emerged, which could not be clearly observed when only looking at single modifications. Indeed, even if a seasonal degree-day factor performed above average in most cases, this particular modification was included in all of the bottom-ranking model structures for model validation against snow water equivalent. Additionally, model structures including an exponential function for snowmelt (M_e) performed above average for all cases and were even included in almost all the top-ranking model structures.

Figure 8Cumulative plots for each of the seven individual modifications to the snow routine of the HBV model as a function of the ranked 64 model structures arising from all the possible combinations of modifications. The modifications are a seasonally variable temperature lapse rate (Γ_s), a linear, sinusoidal, and exponential function for the precipitation phase partition (ΔP_l, ΔP_s, and ΔP_e respectively), different thresholds for precipitation phase and snowmelt (T_P,M), a seasonally variable degree-day factor (C_0,s), and an exponential snowmelt with no refreezing (M_e). Model simulations against snow water equivalent are presented in (a, c), while those against logarithmic stream runoff are presented in (b, d). Model calibration is presented in (a, b) and model validation in (c, d). Model modifications plotted above the 1:1 line (grey dotted line) tend to be included in high-ranking model structures, while those plotted below the 1:1 line tend to be included in low-ranking structures.

Download

Based on the ranked alternative model structures and the modifications contained in each of them, specific model modifications contributed the most to model performance increases, regardless of other model structure modifications. Nevertheless, these dominant modifications impacted the model structure performance in different ways, depending on the modelling scenario (i.e. calibration/validation, objective function). Ideally, any model structure modifications should convey an improved representation of snow water equivalent but also have a positive impact on the simulation of stream runoff (which is the main output of the model), both for model calibration and validation.

To achieve an improved representation of both snow water equivalent and stream runoff, we only took into account those model structures that led to a positive impact for each of the four modelling scenarios (i.e. calibration and validation efforts against both objective functions). This way, we selected and ranked all model alternatives that had a positive median relative model performance value with respect to the HBV model and examined which modifications led to the largest model performance improvement (Fig. 9). All of the selected model structures contained an exponential function for snowmelt (M_e), and none of them included an exponential function for precipitation phase partition (ΔP_e). Most model structures were the result of the combination of three to four individual model structure modifications (seven model structures each). Four model structures contained two model modifications and two model structures contained five modifications. Perhaps most interestingly, two of the model structures included only a single modification: an exponential snowmelt function and a sine function for precipitation phase partition (ΔP_s). Nevertheless, these alternatives had the lowest ranking amongst the selection. Overall, the top-ranking alternatives contained a seasonally varying degree-day factor (C_0,s) and an exponential snowmelt function, while other individual modifications resulted in more considerable model performance variability.

Figure 9Ranked alternative structures of the snow routine of HBV that present positive relative model performance values with respect to HBV for model calibration and validation against snow water equivalent and stream runoff disaggregated by snow routine component variants. The alternatives for each of the considered model components are a linear and seasonally variable degree-day factor (Γ_c and Γ_s respectively); an abrupt, linear, sinusoidal, and exponential precipitation phase partition (ΔP_a, ΔP_l, ΔP_s, and ΔP_e respectively); a common and individualised threshold temperature for precipitation phase partition and snowmelt (T_T and T_P,M respectively); a constant and seasonally variable degree-day factor (C_0,c and C_0,s respectively); and a linear and exponential (with no refreezing) melt function (M_l and M_e respectively). Every row contains one model structure with the selected variant for each of the components highlighted in blue. The median relative model performance for all modelling scenarios is given on the left y axis, while the number of model modifications in each alternative is provided on the right y axis.

Download

While not shown explicitly in this paper, the results obtained from the Monte Carlo sensitivity analysis showed that, even if some of the model structure variants (i.e. T_p,m, ΔP_e, or C_0,s) produced compensating effects on some of the model parameters (e.g. precipitation lapse rate, maximum storage in soil box, threshold for reduction of evaporation, or the shape coefficient), this effect was only observed for a reduced number of catchments. Most parameters showed no compensating effects at all. Overall, parameter values and their sensitivity tended to be reasonably consistent across all the tested model structure variants for most of the catchments in the study.

It is challenging to improve existing rainfall-runoff models and especially those that, like HBV, have successively been tested and applied in many catchments over a range of environmental and geographical conditions (Bergström, 2006). Nevertheless, some of the proposed alternative snow routine model structures that we tested in this study showed a generally positive impact on model performance for simulating both snow water equivalent and stream runoff, albeit to different extents. We found that the most valuable modification to single HBV snow routine components for rainfall-runoff modelling in mountainous catchments in central Europe was the use of an exponential snowmelt function and, to a lesser extent, a seasonally varying degree-day factor. Another modification, using different thresholds for snowfall and snowmelt instead of a single threshold, produced a significant model performance improvement regarding snow water equivalent but did not convey any advantage for simulating stream runoff.

We observed a significant difference in model performance changes between both objective functions when testing the different snow routine model structures. Indeed, in general, the impact was more evident when simulating snow water equivalent than when simulating stream runoff, as the latter is the result of the combined model routines (i.e. snow, soil, groundwater, and routing routines), which partially compensate and mask any modifications made to the snow routine (Clark and Vrugt, 2006). Additionally, some of the modifications that improved the model performance against snow water equivalent, such as the use of an exponential function to define the solid and liquid phases of precipitation, resulted in poorer stream runoff simulations.

Unlike most modifications considered in this study, which are simple conceptualisations of complex processes, the use of an exponential function to describe precipitation phase partition and the use of a seasonally varying temperature lapse rate are both formulations derived from empirical evidence (Magnusson et al., 2014; Rolland, 2002). Nevertheless, as we have previously discussed, neither of these modifications translates into an improvement of model performance for either objective function. This might be because, since models such as HBV are based on simplifications and generalisations of the processes that occur in reality, formulations based on accurate measurements of diverse processes do not align well with the other simplifications made in the model structure and the behaviour at the chosen spatio-temporal resolution (Harder and Pomeroy, 2014; Magnusson et al., 2015).

Other modifications are relatively similar to each other, such as the case of using linear and sine functions to describe precipitation phase partition. Both these formulations require only one additional parameter and perform almost identically: the precipitation partition between rain and snow is exactly the same for both formulations for most of the transition temperature range except for the tails, which are abrupt for the linear case and smooth for the sine one. Provided that the smooth transition is a more accurate description of the physical process which, in addition, avoids the introduction of discontinuities into the objective functions – which might complicate model calibration (Kavetski and Kuczera, 2007), and that both modifications include the same number of parameters and perform nearly identically, the most accurate description should be preferred. Nevertheless, some models, including HBV, continue to use the linear conceptualisation with the argument of simplicity.

Even if we did not observe differences in model performance as a function of different catchment characteristics, snow water equivalent tends to be underestimated in lowland catchments, while we observed no clear pattern for Alpine catchments. Nevertheless, the limited number of high-elevation catchments in the dataset (only four of them have a mean elevation above 2000 m a.s.l.) combined with the generally small size, steep topography, relatively large glacierised areas, scarce vegetation, and exposure to extreme weather conditions (such as strong wind gusts) of these catchments makes it difficult to extract any relevant trends. That being said, in general, the different model structures tended to underestimate snow accumulation and delay the timing of the spring snowmelt season. Similar patterns have also been observed for other mountainous areas of central Europe (Sleziak et al., 2020). We observed differences in model performance among the two geographical domains included in this study.

Among the different snow routine model structures, the modifications to precipitation phase partition penalised the model performance on most Czech catchments for simulating snow water equivalent while having the opposite effect for Swiss catchments. The Czech catchments have a narrower elevation range compared to the Swiss catchments, in addition to an earlier and shorter snowmelt period. These characteristics may favour the simplification of an abrupt transition between rain and snow, while using gradual transitions between rain and snow might favour the more extended melt season and larger elevation ranges of the Swiss catchments. Another factor that may impact the results is the significant differences in model driving and validation data availability for each of the geographical domains (Günther et al., 2019; Meeks et al., 2017). Indeed, while in Czechia there were a limited number of meteorological stations providing temperature, precipitation, and snow water equivalent data, the Swiss catchments benefited from distributed data for the different catchment elevations, allowing for more accurate calibration of snow-related parameters. The difference in resolution between the Swiss and Czech input data might affect the obtained results, where each has its strengths and weaknesses. For instance, high-resolution data can become highly uncertain for individual grid cells, while observational data may be affected by measurement errors and representativeness issues. Even so, the model performance variability of the different snow routine structures relative to the default HBV model should be similar for both cases. Overall, even with the large differences in hydrological regime, catchment morphology, and data availability between the two geographical domains, the impact of the different snow routine model structures on model performance was in general comparable among them.

Regarding model complexity and uncertainty, we found that increasing the degree of detail, and thus the number of parameters of the model, generally translated into a broader range of model performance values, indicating that the uncertainty related to the model structure increased as well. This is a well-known problem of conceptual rainfall-runoff models and the focus of many studies (Essery et al., 2013; Strasser et al., 2002). Additionally, we found that, for most cases, the median model performance increase with an increasing degree of detail was not significant with respect to the performance range. This means that a more detailed model does not necessarily translate into better model performance, which is consistent with previous studies (Orth et al., 2015). This fact highlights the importance of carefully choosing the degree of detail of the model based on the desired objectives and available data (Hock, 2003; Magnusson et al., 2015). Another important aspect is the uncertainty and robustness of the model's parameters. In models, such as the HBV model, model parameters can compensate for each other, which makes the interpretation of any model structure modifications rather challenging (Clark and Vrugt, 2006). Nevertheless, a Monte Carlo sensitivity analysis of the HBV model parameters showed that parameter values and sensitivity were consistent for most catchments and model structures in this study.

Even if an increased degree of detail in the process description is not desirable by itself, as it can lead to overparameterisation and equifinality issues, it can also improve the performance of rainfall-runoff models if it responds to specific needs and available data, among others. Indeed, 22 of the 63 model structure alternatives that we tested in this study (all of them conveying an increase in model detail and therefore an increase in the number of parameters) convey and increase model performance with respect to HBV for both model calibration and simulation against both objective functions. Furthermore, out of these 22 alternatives, only 2 of them consist of a single model structure modification, while most have three or four modifications. Nevertheless, all of these alternatives share common traits, such as using an exponential snowmelt function with no refreezing. Almost all model structures that do not have this particular modification perform worse than HBV in at least one simulation scenario.

It is reasonable to state that, while the increased degree of detail arising from the interplay among the different model structure modifications plays a role in improving the model performance, this is mainly the result of a few dominant modifications. This way, the use of an exponential snowmelt function is the most valuable single modification, with a median performance increase of 0.002 for all simulations (with individual performance increases over 0.1). However, when we combine it with a seasonal degree-day factor, we achieve a median performance increase of 0.008, almost the highest performance increase amongst all model alternatives. Adding further detail to the model does not convey significant improvements for this model structure. Consequently, if we were to implement any modifications to the model, they would be to substitute the linear snowmelt and refreezing conceptualisation with an exponential snowmelt function and replace the constant degree-day factor with a seasonally varying one, in that order.

Finally, it is important to mention that these results are only valid for the selected study areas and cannot be extrapolated to all the different alpine and snow-covered regions around the world as the different processes involved in different geologic, geographic, climatological, and hydrological settings are likely to favour different formulations of the snow processes.

Meteorological and hydrological data to calibrate the HBV model were obtained from the Swiss Federal Office of Meteorology and Climatology, the Swiss Federal Office for the Environment, and the Czech Hydrometeorological Institute. The HBV model outputs are available from the first author upon request.

The supplement related to this article is available online at: https://doi.org/10.5194/hess-24-4441-2020-supplement.

JS initiated the study. MGL developed the methodology and performed all analyses. MGL, NG, and MJ prepared the input meteorological and hydrological data used to calibrate the HBV model. MJPV performed all the necessary modifications to the source code of the HBV model. MGL prepared the manuscript with contributions from all the co-authors.

The authors declare that they have no conflict of interest.

This paper was edited by Elena Toth and reviewed by Giacomo Bertoldi, Juraj Parajka, María José Polo, and Thomas Skaugen.

This research has been supported by the Swiss Federal Office for the Environment (FOEN) and the Czech Science Foundation (project no. GA18-06217Y).