## Abstract

Floods are one of the major climate-related hazards and cause causalities and substantial damage. Accurate and timely flood forecasting and design flood estimation are important to protect lives and property. The Distance Distribution Dynamic (DDD) is a parsimonious rainfall-runoff model which is being used for flood forecasting at the Norwegian flood forecasting service. The model, like many other models, underestimates floods in many cases. To improve the flood peak prediction, we propose a dynamic river network method into the model. The method is applied for 15 catchments in Norway and tested on 91 flood peaks. The performance of DDD in terms of KGE and BIAS is identical with and without dynamic river network, but the relative error (RE) and mean absolute relative error (MARE) of the simulated flood peaks are improved significantly with the method. The 0.75 and 0.25 quantiles of the RE are reduced from 41% to 23% and from 22% to 1%, respectively. The MARE is reduced from 32.9% to 15.7%. The study results also show that the critical support area is smaller in steep and bare mountain catchments than flat and forested catchments.

## INTRODUCTION

Floods are one of the major climate-related hazards and cause causalities and substantial damage on a global scale every year (Hirabayashi *et al.* 2013; Blaikie *et al.* 2014; Winsemius *et al.* 2015). Floods usually cause damage to agricultural land, infrastructure and buildings (Razi *et al.* 2010). Flood peak is one of the most important variables to be estimated as its magnitude and duration are responsible for the damage (Formetta *et al.* 2017; Gao *et al.* 2017). An accurate estimate of flood peak is a critical requirement for proposing appropriate flood damage mitigation measures in order to reduce social and economic costs (Plate 2009).

The common hydrological tools for flood risk management are flood forecasting models and models used to estimate design floods (Plate 2009). The design flood, where the magnitude of the flood is associated with a return period and hence a level of risk, is important in the planning, design and operation of hydraulic structures and for protection of human life and property (Rahman *et al.* 1998; Reis & Stedinger 2005; Smithers 2012). Methods to estimate design floods are generally classified into three: (a) statistical flood frequency analysis; (b) event-based simulation; and (c) derived flood frequency simulation (Filipova *et al.* 2018).

Derived flood frequency analysis, using continuous rainfall-runoff models, is increasing in use for design flood estimation (Cameron *et al.* 2000; Calver & Lamb 1995; Boughton & Droop 2003; Eschenbach *et al.* 2008). A rainfall-runoff model can be used to simulate several flow values under different conditions for extending and enhancing the observed flow record (Filipova *et al.* 2018). A stochastic weather generator is used to simulate long synthetic series of rainfall and temperature input data for the continuous simulation method. The long series of flow data derived from the simulation is then used to estimate the required return periods, usually using plotting positions (Camici *et al.* 2011; Haberlandt & Radtke 2014). There is a growing interest in continuous simulation method of flood estimation as an alternative to event-based method, and internationally the trend is to adopt the continuous method (Lamb & Kay 2004; Chetty & Smithers 2005; Pathiraja *et al.* 2012). The main advantages of the continuous simulation models are their ability to represent the antecedent moisture condition in the catchment and their capability to model future land use and climate changes impacts on the flood peaks (Brocca *et al.* 2011; Smithers *et al.* 2013). The other reason for using the continuous simulation approach is that precipitation records are more widely available and tend to have longer periods of records than stream flow data (Blazkova & Beven 2002). Continuous simulation can avoid the base flow estimation problem in the event-based method and avoid any need to associate return period of the flood with specific design precipitation because the frequency analysis of floods can be done directly.

Rainfall-runoff models are simplified representations of a complex physical system and therefore carry a certain amount of uncertainty in their applications (Bourdin *et al.* 2012). The performance of rainfall-runoff models depend on several factors which include the quality of precipitation input data and an appropriate model structure capable of simulating floods (Collier 2007). Therefore, the structure and performance of the rainfall-runoff models should be evaluated and improved for their capability in simulating flood peaks before using them in design flood estimation and flood forecasting.

There are several ways to classify rainfall-runoff models (Singh 1995). Rainfall-runoff models can be classified into lumped and distributed models. Lumped models consider the whole catchment as a single unit with state variables that represent the average of the catchment (Beven 2001b). Distributed models make prediction at distributed locations, i.e., by discretizing the catchment into a number of elements with state variables representing local averages (Singh & Frevert 2006). When a rainfall-runoff model is used for design flood estimation, the model could underestimate the design flood. Thomas (1982) evaluated floods estimated by continuous simulation methods on 50 small streams in Oklahoma, and the result showed that the flood peaks were consistently underestimated. Pathiraja *et al.* (2012) used 45 catchments in the Murray–Darling basin in Australia to estimate design floods using the Australian water balance model. They found that the model underestimates the floods from 5% to 30% depending on how reasonably the antecedent moisture condition is simulated. The forecast of floods requires an accurate understanding of catchment characteristics and a precise determination of catchment's initial conditions before flooding (Rusjan *et al.* 2009).

*et al.*2011). Gupta

*et al.*(1980) pointed out that the Geomorphic Instantaneous Unit Hydrograph (GIUH) is equivalent to the probability density function of travel times, , from any point in the catchment to the outlet. This permits the formulation of hydrologic response through the geomorphologic width function, . The GIUH and concepts represent the dependency of peak flows on the geomorphological properties of a catchment and provide a quantitative prediction of peak flows for engineering application (Rinaldo

*et al.*1991; Rinaldo

*et al.*1995; D'Odorico & Rigon 2003; Rigon

*et al.*2011). The form and extent of the stream network reflect the characteristics of the hillslope (Willgoose

*et al.*1991). The stream reflects the ground water dynamics and is often termed as perennial, intermittent and ephemeral streams (Dingman 1978; Bencala

*et al.*2011). Dynamic expansions and contractions of stream networks play an important role for hydrologic processes since they connect different parts of the catchment to the outlet (Nhim 2012). Stream networks in a catchment expand and contract as the catchment wets and dries, both seasonally and in response to individual precipitation events, and this dynamic of stream networks gives an important information to the pattern and process of runoff generation (Godsey & Kirchner 2014; Ward

*et al.*2018). The mean of the distribution of distances from a point in the catchment to the nearest river reach ( and the drainage density are among the indexes used to describe a stream network. Horton (1945) defined the traditional as the sum of lengths () of all streams in a catchment divided by the catchment area (). The reflects the spatial characteristics that affect the formation of streams and the response time of a catchment for a particular stream network (Wharton 1994; Tucker

*et al.*2001; Di Lazzaro

*et al.*2014; Skaugen & Onof 2014). The mean distance one has to walk from a random location in a hillslope before encountering a stream, is related to the traditional definition of (Horton 1932; Tucker

*et al.*2001; Di Lazzaro

*et al.*2014).

*et al.*1998) as shown in Equation (2). where

*k*and

*n*are positive numbers. If we insert the value of from Equation (2) into Equation (1), we will get a power relationship between and as shown in Equation (3). where .

The Distance Distributions Dynamics (DDD) model is a parsimonious continuous rainfall-runoff model with a small number of calibration parameters recently developed by Skaugen & Onof (2014). Many of the model parameters can be estimated from catchment topography using GIS and recession characteristics. DDD is a semi-distributed model, i.e., lumped in model parameters and distributed input data (precipitation and temperature). The calibration and validation results for 41 small rural unregulated catchments in Norway (area <50 km^{2}) with hourly data showed that the DDD model, in most cases, underestimated flood peaks (Tsegaw *et al.* 2019). In the runoff dynamics of the DDD model, there is a single static river network forming the basis for the dynamics of water routing through the hillslopes and in the river network (Skaugen & Onof 2014). However, studies show that the river network has a dynamical nature, being more dense during high flows than at low or medium flows (Godsey & Kirchner 2014). The primary objective of this study is to investigate whether including a dynamic river network model into the DDD model will improve flood prediction in small rural catchments (area <50 km^{2}). The secondary objective is to improve the understanding of the stream development for different vegetation covers, catchment slopes and climate. The secondary objective helps us to assess whether there is a potential to relate a calibration parameter of the dynamic river network routine with the environmental factors so that there is a possibility for regionalizing the parameter.

## METHODOLOGY

### Study catchments and data

Fifteen gauged small rural catchments, which show significant underestimation of peak floods during the calibration and validation of 41 small catchments located in Norway (Tsegaw *et al.* 2019), are used in this study for testing the dynamic river network method. The catchments are selected from the Norwegian Water Resources and Energy Directorate (NVE) HYDRA II database. Figure 1 shows locations of the study catchments, and Table 1 shows the catchment characteristics. The definition of small catchment follows that of Fleig & Wilson (2013) with an upper area limit of 50 km^{2}. We selected where the DDD model had a known history of underestimating floods so that the dynamic river network model could be tested and evaluated.

S.no | Cat_ID | Area (km^{2}) | Lake (%) | Marsh land (%) | Forest (%) | Bare mountain (%) | Cultivated land (%) | Urban (%) | Mean elvation (m) | Mean annual precipitation (mm) | Mean annual temprature (°C) | Specific runoff (l/s/km^{2}) | Mean hillslope slope (%) |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|

1 | 6.10 | 7 | 2.7 | 1.6 | 94.3 | 0 | 0 | 0 | 302.9 | 886 | 4.1 | 20.6 | 18.3 |

2 | 12.193 | 50.7 | 1.3 | 3 | 88.4 | 0 | 4.3 | 0.5 | 306.5 | 840 | 3.8 | 17.5 | 15.3 |

3 | 19.107 | 41.5 | 4.8 | 1.4 | 86.4 | 0 | 4.2 | 0.8 | 88.4 | 1,158 | 4.2 | 24.2 | 14.7 |

4 | 26.64 | 9.7 | 7.2 | 2.2 | 38.8 | 46.2 | 1.1 | 0 | 203.9 | 1,688 | 6.7 | 45.8 | 28.3 |

5 | 36.32 | 20.9 | 2.5 | 0.9 | 13.5 | 81.4 | 0 | 0 | 1,039.4 | 2,377 | 0.2 | 105.2 | 34.1 |

6 | 41.8 | 27.4 | 7.7 | 0.4 | 8.8 | 82.2 | 0 | 0 | 836.8 | 2,955 | 2.7 | 126.4 | 37.5 |

7 | 42.2 | 31.1 | 2.8 | 1.4 | 40.7 | 52.1 | 0 | 0 | 573.3 | 2,361 | 4.8 | 108.3 | 40.4 |

8 | 55.4 | 50.6 | 3.7 | 1.3 | 51.8 | 30.7 | 2.8 | 0.2 | 361.3 | 2,593 | 5.4 | 100.4 | 41.9 |

9 | 63.12 | 12.8 | 4.5 | 0.3 | 5.8 | 86.1 | 0.6 | 0 | 886.2 | 2,579 | 1.1 | 94.8 | 34.4 |

10 | 68.2 | 20.9 | 4.3 | 0 | 20.2 | 50.3 | 2.5 | 0 | 402.4 | 2,736 | 5 | 125.2 | 43.6 |

11 | 73.21 | 25.8 | 7.4 | 2.2 | 2.2 | 88 | 0 | 0 | 1,292.5 | 946 | −0.8 | 34.7 | 21.5 |

12 | 73.27 | 30.4 | 9.1 | 0.5 | 0.1 | 89.4 | 0 | 0 | 1,372 | 679 | −2.5 | 33.5 | 14.8 |

13 | 91.2 | 25.8 | 10.1 | 3 | 3.7 | 66.5 | 1.3 | 0 | 261.9 | 2,072 | 6.1 | 63.5 | 29.9 |

14 | 101.1 | 40 | 11.3 | 1.6 | 61.3 | 11.3 | 6 | 0.1 | 232.9 | 1,704 | 5.5 | 54.9 | 23.9 |

15 | 172.8 | 21.2 | 10.7 | 0.5 | 1.4 | 82.5 | 0 | 0 | 659 | 1,465 | 1 | 46 | 17 |

Minimum | 7 | 1.3 | 0 | 0.1 | 0 | 0 | 0 | 88.4 | 679 | −2.5 | 17.5 | 14.7 | |

Maximum | 50.7 | 11.3 | 3 | 94.3 | 89.4 | 6 | 0.8 | 1,372 | 2,955 | 6.7 | 126.4 | 43.6 |

S.no | Cat_ID | Area (km^{2}) | Lake (%) | Marsh land (%) | Forest (%) | Bare mountain (%) | Cultivated land (%) | Urban (%) | Mean elvation (m) | Mean annual precipitation (mm) | Mean annual temprature (°C) | Specific runoff (l/s/km^{2}) | Mean hillslope slope (%) |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|

1 | 6.10 | 7 | 2.7 | 1.6 | 94.3 | 0 | 0 | 0 | 302.9 | 886 | 4.1 | 20.6 | 18.3 |

2 | 12.193 | 50.7 | 1.3 | 3 | 88.4 | 0 | 4.3 | 0.5 | 306.5 | 840 | 3.8 | 17.5 | 15.3 |

3 | 19.107 | 41.5 | 4.8 | 1.4 | 86.4 | 0 | 4.2 | 0.8 | 88.4 | 1,158 | 4.2 | 24.2 | 14.7 |

4 | 26.64 | 9.7 | 7.2 | 2.2 | 38.8 | 46.2 | 1.1 | 0 | 203.9 | 1,688 | 6.7 | 45.8 | 28.3 |

5 | 36.32 | 20.9 | 2.5 | 0.9 | 13.5 | 81.4 | 0 | 0 | 1,039.4 | 2,377 | 0.2 | 105.2 | 34.1 |

6 | 41.8 | 27.4 | 7.7 | 0.4 | 8.8 | 82.2 | 0 | 0 | 836.8 | 2,955 | 2.7 | 126.4 | 37.5 |

7 | 42.2 | 31.1 | 2.8 | 1.4 | 40.7 | 52.1 | 0 | 0 | 573.3 | 2,361 | 4.8 | 108.3 | 40.4 |

8 | 55.4 | 50.6 | 3.7 | 1.3 | 51.8 | 30.7 | 2.8 | 0.2 | 361.3 | 2,593 | 5.4 | 100.4 | 41.9 |

9 | 63.12 | 12.8 | 4.5 | 0.3 | 5.8 | 86.1 | 0.6 | 0 | 886.2 | 2,579 | 1.1 | 94.8 | 34.4 |

10 | 68.2 | 20.9 | 4.3 | 0 | 20.2 | 50.3 | 2.5 | 0 | 402.4 | 2,736 | 5 | 125.2 | 43.6 |

11 | 73.21 | 25.8 | 7.4 | 2.2 | 2.2 | 88 | 0 | 0 | 1,292.5 | 946 | −0.8 | 34.7 | 21.5 |

12 | 73.27 | 30.4 | 9.1 | 0.5 | 0.1 | 89.4 | 0 | 0 | 1,372 | 679 | −2.5 | 33.5 | 14.8 |

13 | 91.2 | 25.8 | 10.1 | 3 | 3.7 | 66.5 | 1.3 | 0 | 261.9 | 2,072 | 6.1 | 63.5 | 29.9 |

14 | 101.1 | 40 | 11.3 | 1.6 | 61.3 | 11.3 | 6 | 0.1 | 232.9 | 1,704 | 5.5 | 54.9 | 23.9 |

15 | 172.8 | 21.2 | 10.7 | 0.5 | 1.4 | 82.5 | 0 | 0 | 659 | 1,465 | 1 | 46 | 17 |

Minimum | 7 | 1.3 | 0 | 0.1 | 0 | 0 | 0 | 88.4 | 679 | −2.5 | 17.5 | 14.7 | |

Maximum | 50.7 | 11.3 | 3 | 94.3 | 89.4 | 6 | 0.8 | 1,372 | 2,955 | 6.7 | 126.4 | 43.6 |

Precipitation, temperature and discharge are the main input data for running and calibrating the DDD model. We used hourly data of precipitation, temperature and discharge. Precipitation and temperature are based on a 1 × 1 km gridded product of the Norwegian Meteorological Institute (http://thredds.met.no/thredds/catalog.html) (Lussana *et al.* 2016). We used a total of 5 years of data for calibration and validation. The DDD model uses distributed precipitation and temperature data as input for the model's 10 elevation zones extracted from the hypsographic curve of a catchment. The elevation of the center of each temperature and precipitation grid cell has been extracted from the 10 × 10 m digital elevation model (DEM) of Norway. Discharge data have been obtained from the Norwegian Water Resources and Energy Directorates (NVE) HYDRA II database. The Norwegian Mapping Authority (www.statkart.no) is the source of the topography, observed river network and land use data.

### The DDD rainfall-runoff model

The DDD model currently runs operationally with daily and three-hourly time steps at the Norwegian flood forecasting service. It has two main modules: the subsurface and the dynamics of runoff. In DDD, the distribution of distances between points in the catchment and their nearest river reach (distance distributions of a hillslope) is the basis for describing the flow dynamics of the hillslope. The distribution of distances between points in the river network and the outlet forms the basis for describing the flow dynamics of the river network. The hillslope and river flow dynamics of DDD is hence described by unit hydrographs (UHs) derived from distance distributions from a GIS and celerity derived from recession analysis (Skaugen & Onof 2014; Skaugen & Mengistu 2016). When the distance distributions are associated with flow celerity of the hillslope and rivers, we obtain the distributions of travel times which constitute the time area concentration curve (Maidment 1993). The derivative of the time area concentration curve gives the instantaneous UH (Bras 1990), which is basically a set of weights distributing the input (precipitation and snowmelt) in time to the outlet.

### Subsurface

*M*(mm), is shared between a saturated zone with volume

*S*(mm) and an unsaturated zone with volume

*D*(mm). If the volume of the saturated zone is high, the unsaturated volume has to be correspondingly small (Skaugen & Onof 2014; Skaugen & Mengistu 2016). The actual water volume present in the unsaturated zone is described as

*Z*(mm). The subsurface state variables are updated after evaluating whether the current soil moisture,

*Z*(

*t*), together with the input of rain and snowmelt,

*G*(

*t*), represent an excess of water over the field capacity,

*R*, which is fixed at 30% (

*R*= 0.3) of

*D*(

*t*) (Skaugen & Onof 2014). If

*G*(

*t*) +

*Z*(

*t*) >

*R**

*D*(

*t*), then the excess water

*X*(

*t*) is added to

*S*(

*t*).

is runoff, and is the actual evapotranspiration which is estimated as a function of potential evapotranspiration and the level of storage. is a degree hour factor which is positive for positive temperature () and zero for negative temperature. is drawn from *Z*. The degree hour approach is a simplification but experiences from Skaugen & Onof (2014) show that the evapotranspiration routine in DDD calculates similar values to the approach used in the well-known rainfall-runoff model HBV (Bergström 1976). A recession analysis of observed runoff from the catchment is used to estimate the catchment-scale fluctuations of storage and the capacity of the subsurface water reservoir (*M*) (see Skaugen & Mengistu 2016).

### Runoff dynamics

The dynamics of runoff in DDD has been derived from the catchment features using a GIS combined with runoff recession analysis. The method for describing the runoff dynamics of a catchment is built on the distance distribution derived from the catchment topography. The distances from the points in the catchment to the nearest river reach are calculated as Euclidean distance for the marsh and soil parts of a hillslope. Previous studies in more than 120 catchments in Norway showed that the exponential distribution describes the hillslope distance (Euclidean distance from the nearest river reach) distribution well, and the normal distribution describes well the distances between points in the river network and outlet of a catchment (Skaugen & Onof 2014). Figure 2 shows the structure of the DDD model. The model is written in the R programming language. All GIS work is done with ArcGIS 10.3 (ESRI 2014), and the recession analysis is done using an R script (R Core Team 2017).

*S*(

*t*), in the catchment (Skaugen & Onof 2014; Skaugen & Mengistu 2016). The celerity associated with the different levels of subsurface storage is estimated by assuming exponential recessions with parameter

*Λ*in the equation , where is the peak discharge immediately before the recession starts. is the slope per of the recession in the log–log space.

The distribution of Λ is modeled using a two-parameter gamma distribution.

*v*, is calculated as a function of Λ using Equation (11). where, is the mean of the distances from points in the catchment (hillslope) to the nearest river. The capacity of the subsurface reservoir

*M*(mm) is divided into five storage levels,

*i*, corresponding to the quantiles of the distribution of Λ under the assumption that the higher the storage, the higher the value of Λ. Each storage level is further assigned a celerity (see Equation (11)), where is the parameter of the UH for the individual storage level

*i,*and estimated such that the runoff from several storage levels will give a UH equal to the exponential UH with a parameter . With the assumption that the recession and its distribution carries information on the distribution of catchment-scale storage, we can consider that the temporal distribution of catchment-scale storage,

*S*(

*t*), is a scaled version to that of Λ.

*S*(

*t*) is calculated using Equation (12), and its distribution is modeled using a two-parameter gamma distribution.

*i*= 1, …, 5). Four storage levels are subsurface level, whereas the fifth one is an overland flow level with unlimited capacity (Skaugen & Onof 2014; Skaugen & Mengistu 2016). The five levels have five-UHs (four for subsurface flow and one for overland flow) and each of them has different temporal scales as they have been assigned different celerities. The UH is modeled as follows: where is the time of input, and is the parameter of the exponential distribution estimated from recession analysis for each level,

*i*.

### Model parameters and calibration

The model parameters are divided into three main groups. The first group are those estimated from observed hydro-meteorological data (Table 2), the second group are those estimated by model calibration against observed discharge (Table 3), and the third group are those estimated from digitized maps using a GIS (Table 4). The snow routine in DDD has two parameters estimated from the spatial distribution of observed precipitation data (Skaugen & Weltzien 2016). The shape parameter (*a*_{0}) and the decorrelation length (*d*) of the gamma distribution of snow and snow water equivalent (SWE) are estimated from a previous calibration for 84 catchments in Norway (Skaugen *et al.* 2015). The calibration of the model is performed using the probability particle swarm optimization (PPSO) algorithm (Lu & Han 2011). The Kling–Gupta efficiency criteria (KGE) have been used as objective function for the calibration (Gupta *et al.* 2009), and we used KGE, the BIAS (ratio of the mean of simulated to observed discharge) and visual inspection of hydrographs to evaluate the performance of the model.

Parameters | Description of the parameter | Method of estimation | Unit |
---|---|---|---|

d | Parameter for spatial distribution of SWE, decorrelation length | From spatial distribution of observed precipitation | Positive real number |

a0 | Parameter for spatial distribution of SWE, shape parameter | From spatial distribution of observed precipitation | Positive real number |

MAD | Long-term mean annual discharge | From long-term observed mean annual flow data | m^{3}/s |

Gshape | Shape parameter of λ | Recession analysis of observed runoff | Positive real number |

Gscale | Scale parameter of λ | Recession analysis of observed runoff | Positive real number |

GshInt | Shape parameter of Λ | Recession analysis of observed runoff | Positive real number |

GscInt | Scale parameter of Λ | Recession analysis of observed runoff | Positive real number |

Parameters | Description of the parameter | Method of estimation | Unit |
---|---|---|---|

d | Parameter for spatial distribution of SWE, decorrelation length | From spatial distribution of observed precipitation | Positive real number |

a0 | Parameter for spatial distribution of SWE, shape parameter | From spatial distribution of observed precipitation | Positive real number |

MAD | Long-term mean annual discharge | From long-term observed mean annual flow data | m^{3}/s |

Gshape | Shape parameter of λ | Recession analysis of observed runoff | Positive real number |

Gscale | Scale parameter of λ | Recession analysis of observed runoff | Positive real number |

GshInt | Shape parameter of Λ | Recession analysis of observed runoff | Positive real number |

GscInt | Scale parameter of Λ | Recession analysis of observed runoff | Positive real number |

Parameters | Description of the parameter | Method of estimation | Unit | Intervals of calibration |
---|---|---|---|---|

pro | Liquid water in snow | Calibration | fraction | 0.03–0.1 |

cx | Degree hour factor for snow melt | Calibration | mm °C^{−1} h^{−1} | 0.05–1.0 |

CFR | Degree hour factor for refreezing | Calibration | mm °C^{−1} h^{−1} | 0.001–0.01 |

Cea | Degree hour factor for evapotranspiration | Calibration | mm °C^{−1} h^{−1} | 0.01–0.1 |

rv | Celerity for river flow | Calibration | m/s | 0.5–1.5 |

Parameters | Description of the parameter | Method of estimation | Unit | Intervals of calibration |
---|---|---|---|---|

pro | Liquid water in snow | Calibration | fraction | 0.03–0.1 |

cx | Degree hour factor for snow melt | Calibration | mm °C^{−1} h^{−1} | 0.05–1.0 |

CFR | Degree hour factor for refreezing | Calibration | mm °C^{−1} h^{−1} | 0.001–0.01 |

Cea | Degree hour factor for evapotranspiration | Calibration | mm °C^{−1} h^{−1} | 0.01–0.1 |

rv | Celerity for river flow | Calibration | m/s | 0.5–1.5 |

Symbol of parameters | Description of the parameter |
---|---|

Area | Catchment area |

maxLbog | Maximum distance of marsh land portion of hillslope |

midLbog | Mean distance of marsh land portion of hillslope |

bogfrac | Areal fraction of marsh land from the total land uses |

zsoil | Areal fraction of DD for soils (what area with distance zero to the river) |

zbog | Areal fraction of distance distribution for marsh land (what area with distance zero to the river) |

midFl | Mean distance (from distance distribution) for river network |

stdFL | Standard deviation of distance (from distance distribution) for river network |

maxFL | Maximum distance (from distance distribution) for river network |

maxDl | Maximum distance (from distance distribution) of non-marsh land (soils) of hill slope |

midDL | Mean distance (from distance distribution) of non-marsh land (soils) of hill slope |

midGl | Mean distance (from distance distribution) for glaciers |

stdGl | Standard deviation of distance (from distance distribution) for glaciers |

maxGl | Maximum distance (from distance distribution) for glaciers |

Hypsographic curve | 11 values describing the quantiles 0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100 |

Symbol of parameters | Description of the parameter |
---|---|

Area | Catchment area |

maxLbog | Maximum distance of marsh land portion of hillslope |

midLbog | Mean distance of marsh land portion of hillslope |

bogfrac | Areal fraction of marsh land from the total land uses |

zsoil | Areal fraction of DD for soils (what area with distance zero to the river) |

zbog | Areal fraction of distance distribution for marsh land (what area with distance zero to the river) |

midFl | Mean distance (from distance distribution) for river network |

stdFL | Standard deviation of distance (from distance distribution) for river network |

maxFL | Maximum distance (from distance distribution) for river network |

maxDl | Maximum distance (from distance distribution) of non-marsh land (soils) of hill slope |

midDL | Mean distance (from distance distribution) of non-marsh land (soils) of hill slope |

midGl | Mean distance (from distance distribution) for glaciers |

stdGl | Standard deviation of distance (from distance distribution) for glaciers |

maxGl | Maximum distance (from distance distribution) for glaciers |

Hypsographic curve | 11 values describing the quantiles 0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100 |

### Dynamic river network routine

We introduce a dynamic river network concept into the DDD model so that the scale of the overland unit hydrograph (OUH) will be dynamic while keeping the four subsurface UHs constant during the simulation period. The methodology, we used in estimating the dynamic OUHs of a hillslope, is similar to that of GIUH and of the width function, i.e., the travel time probability density function of a unit amount of water draining from a catchment. However, the approaches used in estimating the parameters of the distribution are different, i.e., the approach in calculating the celerity and distances of a flow from the points in the hillslope to the nearest river reach. Further, we assumed that the scale of the travel time distribution in a hillslope is dynamic for generating dynamic OUHs while the shape is held constant. In DDD, the dynamic OUHs are turned on and off according to saturation of the subsurface thus giving a dynamic travel time distribution.

The river network indicates where the subsurface water flow becomes surface water flow. The network system governs the dynamics of runoff for conditions where we have no overland flow from the hillslope in that there is a significant (orders of magnitude) difference in water celerity for flow through the soils and flow in the river network (Robinson *et al.* 1995). In case of overland flow, however, we can imagine a dynamic river network (and hence dynamic distance distributions) as a function of overland flow (). We made three assumptions to derive such algorithm.

The mean celerity of the overland flow (

*ν*_{OF}) is constant, i.e., independent of the subsurface saturation and river network.The overland flow unit hydrograph (OUH) is exponential determined from and

*ν*_{OF}.The of a river network is a function of volume of water per area per unit time, i.e., . If we assume a critical flux, of 10 m

^{3}/h is necessary to create a stream, then of 10 mm for an hour over 1000 m^{2}, will provide such a flux, whereas the same flux is obtained for of 100 mm over 100 m^{2}. The two cases have different critical supporting area (), and these cases will provide us with two different river networks where the latter has smaller than the former.

The physical mechanisms underpinning the above three assumptions are:

The variable contributions of saturation excess overland flow of a hillslope develops along the existing river network following the concept of Dunne's overland flow (Dunne 1978).

The critical supporting area ( defines the minimum catchment area from which the generated runoff is sufficient to initiate and maintain river development (Schaefer

*et al.*1979)). The expansion and contraction of the stream network is governed by the amount of saturation excess overland flow.The hillslope travel time probability density function of overland flow is estimated from the distance distributions at any point from the hillslope to the river network and the celerity of flow in the hillslope (D'Odorico & Rigon 2003; Rigon

*et al.*2011).

*ν*

_{OF}. Using assumption (3), we can derive a dynamic after introducing a critical flux as shown in Equation (14), which needs to be determined. where is saturation excess overland flow and is estimated from the DDD model output at each simulation time step. When the subsurface is saturated and there is overland flow (, the dynamic river network subroutine is activated and the corresponding will be calculated in the model using Equation (14).

We need to compute the coefficients *a* and *b* of the general power relation between and of each of the study catchments (see Equation (3)) for computing a dynamic during simulation. For computing *a* and *b*, we have used the following procedure:

The 10 × 10 m DEM of a catchment has been reconditioned to the observed river network using the DEM reconditioning tool from Arc Hydro and a raster flow accumulation map has been prepared using GIS.

We wrote a python script that can loop through several thresholds of flow accumulation () to define stream and create several stream networks. From the distance distributions derived from each stream networks, the is calculated.

A regression curve is fitted to the synthetically derived and of a catchment to estimate and

*b*(Figure 3 shows the fitted curves for six sample catchments). The values of*a*and*b*are unique for each catchment and are listed in Table 5 for all study catchments.

Cat_ID | a | b | R^{2} |
---|---|---|---|

6.1 | 2.17 | 0.37 | 0.99 |

12.193 | 1.16 | 0.42 | 1 |

19.107 | 1.04 | 0.43 | 1 |

26.64 | 1.23 | 0.42 | 1 |

36.32 | 1.18 | 0.42 | 1 |

41.8 | 0.61 | 0.51 | 0.99 |

42.2 | 0.92 | 0.45 | 1 |

55.4 | 1.1 | 0.44 | 1 |

63.12 | 1.05 | 0.45 | 1 |

68.2 | 0.86 | 0.49 | 1 |

73.21 | 0.82 | 0.48 | 1 |

73.27 | 1.19 | 0.44 | 0.99 |

91.2 | 0.68 | 0.51 | 0.99 |

101.1 | 1.1 | 0.45 | 0.98 |

172.8 | 1.08 | 0.44 | 1 |

Cat_ID | a | b | R^{2} |
---|---|---|---|

6.1 | 2.17 | 0.37 | 0.99 |

12.193 | 1.16 | 0.42 | 1 |

19.107 | 1.04 | 0.43 | 1 |

26.64 | 1.23 | 0.42 | 1 |

36.32 | 1.18 | 0.42 | 1 |

41.8 | 0.61 | 0.51 | 0.99 |

42.2 | 0.92 | 0.45 | 1 |

55.4 | 1.1 | 0.44 | 1 |

63.12 | 1.05 | 0.45 | 1 |

68.2 | 0.86 | 0.49 | 1 |

73.21 | 0.82 | 0.48 | 1 |

73.27 | 1.19 | 0.44 | 0.99 |

91.2 | 0.68 | 0.51 | 0.99 |

101.1 | 1.1 | 0.45 | 0.98 |

172.8 | 1.08 | 0.44 | 1 |

After we have obtained the coefficients of the relationship between and , the estimated from Equation (14) will be used to estimate using Equation (3). We estimated from the exponential distribution with parameter , as a distance where 99% of the catchment area is accounted for. From the recession analysis for estimating celerities, we already have an estimate of *ν*_{OF} in the DDD model (Skaugen & Onof 2014). We estimate a dynamic OUH for every time step when overland flow is estimated.

When the calculated using Equation (3) is greater than the of the observed river network, the dynamic river network degenerates to the observed river network. The observed river network is the basis network for all cases where the subsurface capacity is unsaturated, saturated but no overland flow and when there is overland flow but not sufficient to expand the observed (existing) stream network. When the subsurface capacity is saturated and there is sufficient , the observed stream network starts to expand. The extent of expansion is determined by the magnitude of the estimated and .

We have tested the performance of the DDD model with and without the dynamic river network routine. We calibrated and validated the DDD model as described in the model parameter and calibration section, and we implemented the dynamic river routine into the model and calibrated . We have calibrated manually after calibrating automatically DDD parameters without dynamic river network. The procedures we have followed in calibrating are as follows:

The parameter is adjusted by trial and error to fit the observed flood peaks, which had been underestimated by DDD without dynamic river network.

We have visually compared the observed flood hydrographs and flood hydrographs simulated with and without dynamic river network.

While calibrating DDD with (with dynamic river network module), the KGE and BIAS values obtained should not be less than the KGE and BIAS values of DDD without dynamic river network method (earlier calibration result).

Using the visual inspection of observed hydrographs, KGE and BIAS, the which fits the observed flood peaks well is taken as a calibrated value.

*N*is the number of flood peaks (91 in this study). We have also analyzed the quantiles of the distribution of relative errors (RE, Equation (16)) of the flood peaks prediction. where is the observed flood peak and is the predicted flood peak with and without dynamic river network.

### Correlation between with environmental factors

We have done a correlation analysis between the parameters and and environmental factors to improve the understanding on how the dynamic river network develops and to assess the potential for relating to environmental factors. The environmental factors included in the correlation analysis are mean annual precipitation, mean hillslope slope, bare mountain and forest land covers of the study catchments. We have used the Pearson correlation coefficient for the analysis. The derived from observed river network has a spatial variation within a catchment; therefore, we have estimated the mean for each of the study catchments before the correlation analysis using the following two steps. First, the of the observed stream network is estimated using GIS. Second, the mean of the observed river network is calculated using Equation (3). To assess the potential for relating with the environmental factors, correlation analysis between the calibrated and the environmental factors has also been performed. We have also done a stepwise method of multiple linear regression between and the four environmental factors mentioned above to see if there is a possibility for regionalizing .

## RESULTS

### Performance of DDD with and without dynamic river network

The calibration and validation results of the DDD model without a dynamic river network show that the model performs satisfactorily with KGE values between 0.55 and 0.9 and BIAS between 0.75 and 1.25. As stated by Thiemig *et al.* (2013), 0.75 ≤ KGE < 0.9 is considered good, 0.5 ≤ KGE < 0.75 is intermediate and 0.0 ≤ KGE < 0.5 is poor. Seven catchments show good and eight catchments show intermediate performance both during the calibration and validation periods. Even if the KGE performance is satisfactory, the visual inspection of the hydrographs shows that several observed flood peaks are underestimated.

We added the dynamic river network routine into the DDD model and calibrated the critical flux () parameter of the routine manually for the whole simulation period. The KGE and BIAS performance of the model are similar as before, i.e., without dynamic river network for all study catchments except one, where the KGE is slightly lower. However, the inspection of the hydrographs clearly shows that the predication of several underestimated flood peaks has been improved after the addition of the dynamic river network routine. The dynamic OUHs that resulted from the dynamic river network have higher peaks and narrower width during the flooding events, and these OUHs, added with the subsurface UHs, helped in improving the previously underestimated floods. Figure 4 shows the hillslope distance distributions for variable for catchment 73.27. Figure 5 shows the empirical cumulative distance distributions functions as an example for the dynamic distance distribution presented in Figure 4, and Figure 6 shows the four dynamic OUHs which resulted from the corresponding distance distributions functions. Table 6 shows , and for a catchment 12.193 during a flooding event.

OF (mm/h) | A_{c} (m^{2}) | D_{m} (m) estimated using Equation (3) | D_{m} (m) used in deriving OUH |
---|---|---|---|

0.144 | 1,319,444.4 | 429.08 | 301.1 |

2.19 | 86,758 | 137.58 | 137.58 |

2.22 | 85,585.6 | 136.79 | 136.79 |

0.94 | 202,127.7 | 196.26 | 196.26 |

10.9 | 17,431.2 | 70.12 | 70.12 |

1.4 | 135,714.3 | 166.02 | 166.02 |

2.7 | 70,370.4 | 126 | 126 |

1.19 | 159,663.9 | 177.75 | 177.75 |

0.19 | 1,027,027 | 388.44 | 301.1 |

1.27 | 149,606.3 | 172.96 | 172.96 |

0.2 | 950,000 | 375.93 | 301.1 |

1.33 | 142,857.1 | 169.64 | 169.64 |

1.28 | 148,902.8 | 172.61 | 172.61 |

0.17 | 1,117,647.1 | 402.48 | 301.1 |

0.14 | 1,319,444.4 | 431.54 | 301.1 |

0.64 | 296,875 | 230.64 | 230.64 |

0.34 | 558,823.5 | 302 | 301.1 |

0.31 | 612,903.2 | 312.72 | 301.1 |

0.83 | 228,915.7 | 206.79 | 206.79 |

OF (mm/h) | A_{c} (m^{2}) | D_{m} (m) estimated using Equation (3) | D_{m} (m) used in deriving OUH |
---|---|---|---|

0.144 | 1,319,444.4 | 429.08 | 301.1 |

2.19 | 86,758 | 137.58 | 137.58 |

2.22 | 85,585.6 | 136.79 | 136.79 |

0.94 | 202,127.7 | 196.26 | 196.26 |

10.9 | 17,431.2 | 70.12 | 70.12 |

1.4 | 135,714.3 | 166.02 | 166.02 |

2.7 | 70,370.4 | 126 | 126 |

1.19 | 159,663.9 | 177.75 | 177.75 |

0.19 | 1,027,027 | 388.44 | 301.1 |

1.27 | 149,606.3 | 172.96 | 172.96 |

0.2 | 950,000 | 375.93 | 301.1 |

1.33 | 142,857.1 | 169.64 | 169.64 |

1.28 | 148,902.8 | 172.61 | 172.61 |

0.17 | 1,117,647.1 | 402.48 | 301.1 |

0.14 | 1,319,444.4 | 431.54 | 301.1 |

0.64 | 296,875 | 230.64 | 230.64 |

0.34 | 558,823.5 | 302 | 301.1 |

0.31 | 612,903.2 | 312.72 | 301.1 |

0.83 | 228,915.7 | 206.79 | 206.79 |

Italic numbers are rounded to two significant figures.

Figure 7 shows the hydrographs during the flooding periods with and without the dynamic river network routine for six sample study catchments. Table 7 shows the observed floods, simulated floods, KGE and BIAS performance of DDD model with and without dynamic river network routine for five sample study catchments selected randomly.

Cat_ID | Observed flood(s) in m^{3}/s | Simulation period | Simulated value(s) of flood without river dynamics in m^{3}/s | Performance of DDD model without river dynamics in calibration | Simulated value(s) of flood with river dynamics in m^{3}/s | Performance of DDD model with river dynamics | ||
---|---|---|---|---|---|---|---|---|

KGE | BIAS | KGE | BIAS | |||||

12.193 | 29.42 | 2 years | 12.97 | 0.64 | 1.2 | 28.78 | 0.65 | 1.2 |

19.107 | 11.65 and 16.2 | 3 years | 7.64 and 10.3 | 0.8 | 0.93 | 9.42 and 16.1 | 0.81 | 0.94 |

41.8 | 43.96 and 36.13 | 2 years | 28.02 and 26.74 | 0.77 | 0.84 | 36.3 and 40.3 | 0.77 | 0.84 |

73.27 | 18.85 | 3 years | 13.3 | 0.71 | 0.76 | 18.3 | 0.71 | 0.76 |

91.2 | 12.06 | 2 years | 8.34 | 0.71 | 0.8 | 12.04 | 0.71 | 0.8 |

Cat_ID | Observed flood(s) in m^{3}/s | Simulation period | Simulated value(s) of flood without river dynamics in m^{3}/s | Performance of DDD model without river dynamics in calibration | Simulated value(s) of flood with river dynamics in m^{3}/s | Performance of DDD model with river dynamics | ||
---|---|---|---|---|---|---|---|---|

KGE | BIAS | KGE | BIAS | |||||

12.193 | 29.42 | 2 years | 12.97 | 0.64 | 1.2 | 28.78 | 0.65 | 1.2 |

19.107 | 11.65 and 16.2 | 3 years | 7.64 and 10.3 | 0.8 | 0.93 | 9.42 and 16.1 | 0.81 | 0.94 |

41.8 | 43.96 and 36.13 | 2 years | 28.02 and 26.74 | 0.77 | 0.84 | 36.3 and 40.3 | 0.77 | 0.84 |

73.27 | 18.85 | 3 years | 13.3 | 0.71 | 0.76 | 18.3 | 0.71 | 0.76 |

91.2 | 12.06 | 2 years | 8.34 | 0.71 | 0.8 | 12.04 | 0.71 | 0.8 |

Italic numbers are rounded to two significant figures.

The results of the statistical analysis (mean absolute relative errors and quantiles of relative errors), for the 91 observed peak floods from the 15 study catchments, show that the dynamic river network method improved the prediction of peak floods significantly. The 0.75 quantile of the relative errors of the simulated peaks reduced from 41% to 23%, and the 0.25 quantile of the relative errors reduced from 22% to 1%. Figure 8 shows box plots of the relative errors with and without river dynamics. The of the magnitude of the underestimated peak floods is reduced from 32.9% to 15.7%.

### Correlation between with environmental factors

The critical supporting area, , of an observed stream network of a catchment shows spatial variation within the catchment (Figure 9 shows the distributions for five sample catchments); therefore, the mean value of a catchment is used for the correlation analysis. The mean for the observed river networks is correlated with environmental factors, i.e., vegetation cover, topography and climate. The correlation with vegetation cover is stronger than that of topography and climate. The mean has a positive correlation with the forest cover in a catchment, but it has a negative correlation with mean annual precipitation, bare mountain cover and mean hillslope slope of a catchment. Table 8 shows the correlation between the mean and and the environmental factors.

Forest (%) | Bare mountain (%) | Mean hillslope slope (%) | Mean annual precipitation (mm) | |
---|---|---|---|---|

Mean A_{c} (m^{2}) | 0.63 | −0.63 | −0.56 | −0.49 |

F_{c} (m^{3}/h) | 0.46 | −0.5 | −0.29 | −0.18 |

Forest (%) | Bare mountain (%) | Mean hillslope slope (%) | Mean annual precipitation (mm) | |
---|---|---|---|---|

Mean A_{c} (m^{2}) | 0.63 | −0.63 | −0.56 | −0.49 |

F_{c} (m^{3}/h) | 0.46 | −0.5 | −0.29 | −0.18 |

The calibrated critical flux, , of the dynamic river network routine is correlated with the environmental factors. The correlation between and the vegetation cover is stronger than the correlation between and topography and mean annual precipitation. shows positive correlation with forest and negative correlation with bare mountain, mean annual precipitation and mean hillslope slope of a catchment. Table 9 shows the environmental factors used in the correlation analysis, the and mean of observed river network and the calibrated of the dynamic river network routine of the study catchments.

Cat_ID | Environmental factors | Observed river network | Dynamic river network | ||||
---|---|---|---|---|---|---|---|

Forest (%) | Bare mountain (%) | Mean annual precipitation (mm) | Mean hillslope slope (%) | D_{m} (m) | Mean A_{c} (m^{2}) | F_{c} (m^{3}/h) | |

6.1 | 94.3 | 0 | 886 | 18.3 | 149.2 | 92,406 | 15 |

12.193 | 88.4 | 0 | 840 | 15.3 | 301.1 | 560,042 | 190 |

19.107 | 86.4 | 0 | 1,158 | 14.7 | 336.9 | 689,123 | 370 |

26.64 | 38.8 | 46.2 | 1,688 | 28.3 | 181.6 | 146,146 | 5 |

36.32 | 13.5 | 81.4 | 2,377 | 34.1 | 168.6 | 135,213 | 120 |

41.8 | 8.8 | 82.2 | 2,955 | 37.5 | 157.2 | 53,418 | 15 |

42.2 | 40.7 | 52.1 | 2,361 | 40.4 | 175.9 | 117,477 | 30 |

55.4 | 51.8 | 30.7 | 2,593 | 41.9 | 155.1 | 76,715 | 15 |

63.12 | 5.8 | 86.1 | 2,579 | 34.4 | 181.2 | 93,596 | 80 |

68.2 | 20.2 | 50.3 | 2,736 | 43.6 | 211.3 | 75,573 | 150 |

73.21 | 2.2 | 88 | 946 | 21.5 | 298.4 | 216,464 | 15 |

73.27 | 0.1 | 89.4 | 679 | 14.8 | 189.4 | 100,964 | 60 |

91.2 | 3.7 | 66.5 | 2,072 | 29.9 | 283.8 | 137,481 | 90 |

101.1 | 61.3 | 11.3 | 1,704 | 23.9 | 334.0 | 328,354 | 150 |

172.8 | 1.4 | 82.5 | 1,465 | 17 | 168.1 | 95,972 | 10 |

Cat_ID | Environmental factors | Observed river network | Dynamic river network | ||||
---|---|---|---|---|---|---|---|

Forest (%) | Bare mountain (%) | Mean annual precipitation (mm) | Mean hillslope slope (%) | D_{m} (m) | Mean A_{c} (m^{2}) | F_{c} (m^{3}/h) | |

6.1 | 94.3 | 0 | 886 | 18.3 | 149.2 | 92,406 | 15 |

12.193 | 88.4 | 0 | 840 | 15.3 | 301.1 | 560,042 | 190 |

19.107 | 86.4 | 0 | 1,158 | 14.7 | 336.9 | 689,123 | 370 |

26.64 | 38.8 | 46.2 | 1,688 | 28.3 | 181.6 | 146,146 | 5 |

36.32 | 13.5 | 81.4 | 2,377 | 34.1 | 168.6 | 135,213 | 120 |

41.8 | 8.8 | 82.2 | 2,955 | 37.5 | 157.2 | 53,418 | 15 |

42.2 | 40.7 | 52.1 | 2,361 | 40.4 | 175.9 | 117,477 | 30 |

55.4 | 51.8 | 30.7 | 2,593 | 41.9 | 155.1 | 76,715 | 15 |

63.12 | 5.8 | 86.1 | 2,579 | 34.4 | 181.2 | 93,596 | 80 |

68.2 | 20.2 | 50.3 | 2,736 | 43.6 | 211.3 | 75,573 | 150 |

73.21 | 2.2 | 88 | 946 | 21.5 | 298.4 | 216,464 | 15 |

73.27 | 0.1 | 89.4 | 679 | 14.8 | 189.4 | 100,964 | 60 |

91.2 | 3.7 | 66.5 | 2,072 | 29.9 | 283.8 | 137,481 | 90 |

101.1 | 61.3 | 11.3 | 1,704 | 23.9 | 334.0 | 328,354 | 150 |

172.8 | 1.4 | 82.5 | 1,465 | 17 | 168.1 | 95,972 | 10 |

## DISCUSSION

### Dynamic river networks

Dynamic river networks and hence dynamic OUHs are introduced and implemented in the DDD rainfall-runoff model to improve the simulation of floods. The dynamic river network method expands the observed river networks during events. The expansion means that the required to initiate and maintain a stream decrease. Smaller results in smaller (see Table 5). The smaller value indicates shorter travel times from points in the catchment to the nearest river reach. The shorter travel time distribution generates OUHs with a higher peak and shorter scale for the hillslopes (Figure 6). The dynamic OUHs are superpositioned with the other four subsurface UHs of DDD to give a single dynamic UH of a catchment during flooding events. The results of the method are supported by the previous study of D'Odorico & Rigon (2003) who found that shorter hillslope distances result in shorter travel times and hence higher flood peaks. The smaller results, obtained during the flooding events using the dynamic river network method, are also supported by the study of Humbert (1990) who found that a good correlation exists between the runoff coefficient of flooding events and , and hence the , for 45 French catchments. Lazzaro *et al.* (2015) also found that the variability of runoff due to higher (lower ) creates a faster concentration of flow that implies shorter travel times and higher peak floods. The results in this study are also supported by the results of Lee *et al.* (2008) who found that a UH of a catchment is dynamic during different precipitation intensities, i.e., the higher the precipitation intensity, the higher the peak and shorter the temporal scale of the UHs. The results of this study also show that a dynamic river network method could be a solution for rainfall-runoff models which face challenges in predicting flood peaks through continuous simulation. Improving the prediction of peak floods in a continuous simulation is very important because the hydrograph consisting of this peak flow is mainly responsible for the damage caused by floods. Therefore, a dynamic river network is a method to be conceptualized and included as one routine in continuous rainfall-runoff models which underestimate predictions of floods.

We analyzed statistically 91 underestimated flood peaks to evaluate the performance of the dynamic river network. The and quantiles of of the prediction with and without dynamic river network show that the overall performance of the method has improved the prediction of the peaks satisfactorily. The dynamic river network overestimated 17 of the 91 flood peaks and still underestimates the remaining 74 floods but with a significant improvement in the prediction of flood peaks compared to the results obtained without a dynamic river network.

A single calibrated critical flux, improves the prediction of several underestimated floods significantly, but it also overestimates a few flood peaks (Figure 7). Reasons for overestimation could be that a single calibrated could not represent the different precipitation patterns, overland flow patterns and initial conditions prior to flooding events. Another reason for such overestimation could be that the dynamic estimated using Equation (3) is not always an accurate representation of the reality during a flood event. The limitations of manual calibration include subjectivity and time consuming. Manual calibration methods are subjective in the evaluation of model fit and the final choice of optimal parameters. Ndiritu (2009) pointed out that manual calibration may be more prone to obtaining suboptimal parameter sets than automatic calibration. Studies also show that manual calibration is more subjective than automatic calibration because it largely depends on visual hydrograph inspection and the personal judgment of the hydrologist. Substantial amount of time is also required to adjust so that the observed and simulated flood peaks agreed well. A separate automatic calibration of (which is not included in this study) after defining and writing appropriate objective function could improve the limitations of the manual calibration. However, since we have only one manually calibrated parameter and we have enough experience of using DDD model, the manual calibration result of could be very good. In addition, the results of our study show that the manually calibrated resulted in a significant improvement in predicting flood peaks using dynamic river network method.

The hydrograph in Figure 7 shows that the two floods of the catchment 41.8 could not be estimated well using a single calibrated of magnitude 15 m^{3}/h even if the overall prediction of the flood peaks is improved. When we look at the flood hydrographs, the 27 November flood of 43.96 m^{3}/s (at 8 A.M), had been preceded by a 1-day precipitation of 68.9 mm (from 26 November 2011 09:00 to 27 November 2011 08:00). The precipitation was again preceded by a 3-day precipitation of magnitude 86.1 mm, i.e., the 4-day precipitation preceded the flooding event was 154.9 mm. When is fitted to this single flooding event, 5 m^{3}/h is required. After 1 month, another heavy precipitation event happened (56.9 mm/day) and the event was preceded by 82.2 mm of 3-day precipitation, i.e., the 4-day precipitation preceded the flooding event was 139.7 mm. The magnitude of the flood was 36.13 m^{3}/s. When is fitted to this single flooding event, 25 m^{3}/h is required. The variation in the fitted values of for different flooding events in a catchment shows that we could have an overestimation of flooding events when we use a single calibrated for the whole simulation period as a representative for a catchment. Accordingly, a single of 15 m^{3}/h for the catchment 41.8 has overestimated the December 2011 flood, i.e., 36.13 m^{3}/s, but it has improved the overall prediction of the flood peaks in the catchment.

The spatial variability of during flooding events, which is not considered in this study, could also be another factor for the overestimation of floods using a single calibration value of for a catchment. We have derived the coefficients of Equation (3) (*a* and *b*) with the assumption of constant using the DEM, which considers only the topography of a catchment. However, the of observed stream networks clearly shows that there is a spatial variability of within a catchment. Geological and land use factors play significant roles in initiating and maintaining a stream, and these factors control the spatial variability of in a catchment (Montgomery & Dietrich 1989; Ogden *et al.* 2011; Sjöberg 2016; Ward *et al.* 2018) in addition to the topography. The correlation results between and vegetation cover, which is done in this study and explained in the next section, also confirm that the land use affects the spatial variability of . Figure 7 shows box plots of the spatial distributions of of observed stream network for five sample catchments. For flooding events preceded by short duration and higher values of , Equation (3) gives very low values of . The very low gives OUHs of sharp peak and short scales which overestimated the floods. However, if we had calculated the actual using the spatial variability of , we could have found higher values of than the value calculated using Equation (3) and the overestimation could have been avoided.

For estimating the parameters of travel time distributions of overland flow of a hillslope, we followed the original approach used in DDD (e.g. the distance from any point in the catchment to the nearest river network is calculated using the Euclidean distance and the celerity is determined from recession analysis). The GIUH and of width function estimates the distribution of travel times at the outlet of a catchment combining the hillslope and river network travel times using the steepest descent path from any point in the catchment to the outlet and the shape and scale parameters of the travel time distribution could change with the extent of hillslope saturation. Therefore, further investigation, i.e., comparison assessment is required before concluding as one method is better than the other.

### Correlation of and with environmental factors

Environmental factors such as vegetation cover, topography and climate, affect and hence . Land use (e.g. vegetation cover) affects the hydrology and can affect subsurface as well as overland flow which in turn can cause changes in the stream network, i.e., (Jonathan & Dennis 2001). The correlation results show that the denser the vegetation cover in a catchment the higher the is required for initiating and maintaining a stream and vice versa, i.e., positive correlation with forest and negative correlation with bare mountain. The correlation result confirms the findings that a decrease in vegetation causes a decrease in surface resistance and critical shear stress, which result in an increase of drainage density (a reduction in ), i.e., streams form easier in less vegetated catchments (Willgoose *et al.* 1991; Prosser & Dietrich 1995; Magnuson *et al.* 1997; Tucker & Slingerland 1997). Field observations also show that higher and hence higher is generally observed in denser vegetation cover (Morisawa 1985).

The steepness of a catchment is one of the topographical factors controlling and hence . In this study, we used the mean hillslope slope of a catchment and found that a catchment hillslope slope has negative correlation with , i.e., the higher the steepness, the lower the required to initiate and maintain stream. This finding is supported by Montgomery & Dietrich (1989) who found that stream initiation on steep slopes shows a negative relationship between valley gradient at the stream head and , i.e., the higher the stream head slope the lower the (lower ).

The positive correlation between and vegetation cover in a catchment shows that the denser the vegetation covers, the higher the . shows negative correlation with bare mountain, mean hillslope slope and mean annual precipitation of a catchment. The and hence depend on several factors, which include geology, precipitation, vegetation, morphology, soils and land uses, and one factor may be more important than another (see Table 8). Therefore, a more detailed investigation supported by field work (e.g. mapping of the slope, geology, vegetation cover and soil of a catchment at the head of first-order streams of observed river networks and mapping of the pattern of expansion of first-order streams during flooding events) should be carried out to assess how the combination of these factors control and hence .

We have done a simple multiple linear regression analysis using the four environmental factors as predictors, i.e., forest, bare mountain, slope and mean annual precipitation, to estimate the response variable . The result shows that only bare mountain is contributing significantly in estimating with a significance level of 0.1, and the coefficient of determination of the regression () is 0.3. The objective of the regression analysis, done in this study, is to assess a preliminary method for regionalization that can predict for ungauged catchments from environmental factors and to lay a foundation for further studies.

A dynamic river network method could be implemented in rainfall-runoff models as shown for DDD for prediction of floods for catchments with a wide range of topography and land uses (Table 1). In this study, the effect of steep hills is reflected in the dynamic river networks as the steepness of a catchment is one of the factors that govern the initiation of streams. As shown in Table 8, the mean hillslope slope has a negative correlation with , i.e., we need a smaller to initiate and maintain streams in steep topography than in a flat topography. Table 9 also shows that the mean of an observed stream network decreases as the mean hillslope slope of a catchment increases. The fundamental theory behind the method is the expansion of river networks during flooding events, i.e., whether the critical flux, , which is required to initiate and maintain a stream, is satisfied or not. The magnitude of depends on the magnitude of saturation excess overland flow, , and the critical support area, . The study results also show that the critical support area, required to initiate and maintain a stream, is smaller in steep and bare mountain catchments than flat and forested catchments. Therefore, the method could be applicable at different catchments with different characteristics.

## CONCLUSIONS

The dynamic river network method, introduced in Distance Distribution Dynamics (DDD) rainfall-runoff model, can improve the prediction of flood peaks in continuous simulation satisfactorily. The performance of the DDD model is the same with and without dynamic river network in terms of KGE and BIAS. The statistical analysis on 91 flood peaks, underestimated by DDD without dynamic river network method, shows that the MARE of the prediction reduced from 32.9% to 15.7% using the dynamic river network method. With a dynamic river network method, the 0.75 quantile of the relative errors reduced from 41% to 23%, and the 0.25 quantile of the relative errors reduced from 22% to 1%. The visual inspection of the hydrographs also shows an improvement in the prediction of flood peaks for several flooding events. Therefore, we recommend the use of a dynamic river network method in the prediction of floods. The next step in the development of the method is to investigate the applicability of the method from gauged to ungauged catchments and find a way to address the limitations identified in this study.

The critical flux, the calibration parameter introduced in the method, has been formulated as the product of critical supporting area and the saturated excess overland flow . shows stronger correlation with vegetation cover than topographical and climate factors. The parameter shows positive correlation with forest cover of catchments, and negative correlation with bare mountain, mean hillslope slope and mean annual precipitation. The simple multiple linear regression, using the four environmental factors as predictors and as a response variable, shows that there is a potential to estimate from environmental factors and regionalize it for using the method without calibration. The value of the calibrated could be different for the same catchment of different flood magnitudes, and it could be different for the same type of vegetation cover for different catchments. This difference shows that could depend on several environmental factors and further investigations should be carried out.

## ACKNOWLEDGEMENTS

The authors would like to acknowledge Cristian Lussana of the Norwegian Meteorological Institute for providing information on how to access and process the 1 × 1 km spatial and 1-h temporal resolution gridded precipitation and temperature data for Norway. Finally, the authors gratefully acknowledge the financial support by the Research Council of Norway and several partners through the Centre for Research-based Innovation ‘Klima 2050’ (see www.klima2050.no).

## REFERENCES

*Geographical Research Letters*

**9**

*Earth Surface Processes and Landforms*

**26**, 115–134