Advances in Water Resources Stream-scale ﬂow experiment reveals large inﬂuence of understory growth on vegetation roughness

Vegetation is a key source of ﬂow resistance in natural channels and ﬂoodplains. It is therefore important to ac- curately model the ﬂow resistance to inform decision makers and managers. However, it is challenging to predict the resistance of real vegetation, because vegetation models are based on relatively small-scale lab experiments with mostly artiﬁcial vegetation. Experimental tests of real vegetation under ﬁeld conditions are scarce. The purpose of this study is to measure the ﬂow resistance of a submerged willow patch, where small herbaceous vegetation was allowed to grow in between the willow stems to simulate ﬁeld conditions. Detailed ﬂow velocity measurements were performed during an full scale experiment of ﬂow around a submerged patch of willows. The parameter values of the willow vegetation model, as well as the friction coeﬃcients of the vegetated banks and unvegetated channel bed, were computed simultaneously using Bayesian inference using a 2D hydrodynamic model. Results show that the presence of understory growth greatly aﬀects ﬂow patterns and the value of the eﬀective vegetation density parameter. Measured ﬂow velocities in the patch with understory growth were very low, and the patch has relatively high deﬂection. After removal of this undergrowth, ﬂow velocities in the patch increased and deﬂection of the vegetation canopy decreased. We show that estimating vegetation density using an often-used rigid cylinder estimator based on vegetation sampling, underestimated the eﬀective value by more than an order of magnitude. We argue that proposed extensions to existing vegetation models, which can take into account understory growth and reconﬁguration, could be tested under ﬁeld conditions using the approach followed in this paper.


Introduction
Vegetation in rivers and streams is one of the largest sources of flow resistance ( Luhar and Nepf, 2013 ) and an important source of uncertainty in hydrodynamic models used for flood management ( Warmink et al., 2013 ) and river engineering applications ( Berends et al., 2018 ). The presence of vegetation affects the morphological evolution of rivers ( van Oorschot et al., 2015 ), biodiversity ( Straatsma et al., 2017 ), and water quality ( Dosskey et al., 2010 ). Therefore, a good representation of vegetation in computer models is important.
Advances in remote sensing and classification algorithms enable detailed maps of the spatial distribution of vegetation ( Geerling et al., 2009;Forzieri et al., 2011;Hodges, 2015 ). The local effect of vegetation on flow is commonly modelled through a friction term in the momentum equation. This term, also known as roughness or flow resistance, often represents various ways of energy losses that are not explicitly modelled. We can broadly distinguish two ways to determine vegetative friction in field situations.
The first approach to determine vegetative friction is based on lookup tables or vegetation models. The table of typical (Manning-type) friction values by Chow (1959) are still used to characterise the roughness of streams based on images or descriptions of the vegetation (e.g. "light brush and trees, in winter "). Look-up tables are a convenient way to couple remote-sensing to hydrodynamic modelling (e.g., Forzieri et al., 2011 ). However, despite the widespread use of the Manning coefficient as a lumped parameter to characterise describe roughness, the implicit with the adopted coordinate system. Right: schematic cross-section at the location of the first patch (crosssection D3).
does not hold for vegetation ( Naden et al., 2006;Ferguson, 2007 ). For this reason, various (semi-) empirical models have been developed that compute the contribution for vegetation on flow resistance based on vegetation parameters, such as stem count and plant morphology ( Klopstra et al., 1997;Baptist et al., 2007;Huthoff et al., 2007;Yang and Choi, 2010;Li et al., 2015 ). It has been shown that such models perform well against a large body of data collected from laboratory experiments ( Vargas-Luna et al., 2015 ).
Further studies have been carried out to better capture real-world vegetation dynamics such as reconfiguration ( Järvelä, 2004;Dijkstra and Uittenbogaard, 2010;Verschoren et al., 2016 ) and contribution of foliage ( Bal et al., 2011;Västilä and Järvelä, 2014 ). In practice, the vegetation parameters necessary to use these formulas should be either directly measured, or taken from a look-up table. Successful application of this approach to field conditions depends on the validity of the vegetation model and the availability of data on vegetation.
The second approach to determine vegetative friction is through inverse modelling, which involves fitting the model (in practice usually a subset of model parameters) to measurements. It is called 'inverse' modelling, because the quantity of interest is not model output, but model input (such as parameters). In practice, inverse modelling is used to find the values for look-up tables, because the lumped roughness can be straightforwardly computed by inverting the Manning equation, under the assumption of steady uniform flow, a logarithmic flow profile and a known energy slope ( Marjoribanks et al., 2014 ). In the case of non-uniform flow, the inverted Bélanger equation can be used instead ( Errico et al., 2018 ).
Roughness values estimated through inverse modelling are often used as the 'measured' roughness to provide the experimental basis for vegetation models. However, estimating the friction in heterogeneous real-world situations based on the energy slope (e.g., a single vegetation patch which does not cover the entire channel), is complicated. Estimating the friction of multiple sources, based on a single observation (the energy slope), may result in non-unique solutions to parameter values ( Beven, 2006 ). This problem is analogous to that of underdetermination in regression problems, which allows an infinite number of acceptable combinations of unknown friction factors if the parameters exceed the available observations. This severely limits our ability to estimate the friction values, and validate vegetation models, in real-world situations.
In cases where non-uniqueness is an issue, probabilistic parameter estimation is preferred to deterministic optimisation ( Matott et al., 2009;Guillaume et al., 2019 ). Formally, this is achieved through Bayesian inference, which produces probability distributions of parameter values given the (subjective) likelihood that the model, given these values, confirms experimental observation. A well known Bayesian inference methodology is GLUE (generalised likelihood uncertainty estimation), originally developed for non-uniqueness problems in hydrology ( Beven and Binley, 1992;Fonseca et al., 2014;Mayotte et al., 2017 ). However, previous work has shown that estimating more than one friction source based on water levels can result in unidentifiable parameter distributions, which are characterised by wide and unconstrained shapes ( Werner et al., 2005 ). To increase the identifiability of parameter distributions other observational data can be used, such as inundation patterns ( Pappenberger et al., 2005 ). Hypothetically, detailed velocity measurements around vegetation patches can serve a similar function to es-timate flow resistance using GLUE. However, the literature on flow data around real-world vegetation patches are scarce ( Naden et al., 2006;Marjoribanks et al., 2017 ), and comparisons between laboratory and field in general are rare ( Huthoff et al., 2013;Groom and Friedrich, 2018 ).
Our aim is to investigate the flow resistance of vegetation under natural conditions. Specifically, we are interested in the effect of secondary vegetation growing in between and under the branches of the dominant vegetation species. Such secondary vegetation is hard to detect remotely from satellite or drone imagery, and may therefore lead to an underestimation of biomass in large-scale river models that rely on ecotope maps to estimate floodplain friction ( Straatsma and Huthoff, 2011;Forzieri et al., 2011;Warmink et al., 2013;van Oorschot et al., 2015;Berends et al., 2019 ). It is expected that understory growth increases the flow resistance, but it is unknown by how much. In this study, we perform a stream-scale experiment with real vegetation where secondary vegetation was allowed to develop naturally. This physical experiment is coupled to a digital twin numerical model to compute friction parameters by Bayesian inference.

Physical experiment setup
Large-scale experiments were performed at the site of the KICT-REC (Korea Institute of Civil Engineering and Building Technology River Experiment Center), which is located in the city of Andong, South Korea. This facility is designed for full scale physical experiments and consists of three separate channels of various slope and sinuosity. Large capacity pumps can generate the maximum flow rate up to 10 m 3 s −1 . The length of each channel is approximately 600 m.
The experiments for vegetated flow were performed in the downstream section of the "A1 " channel, which has a trapezoidal cross-section with a bottom width of 3 m, top width of 11 m, bed slope of 0.001 m/m and bank-side slope of roughly 1:1.5 (V:H), as shown in Fig. 1 . The bankside slopes of the channel are vegetated with grasses and small annuals native to the region. The channel bed consists of sandy materials with a mean particle size of 0.8 mm.
Seven alternating willow patches, each with a length of 4 m and width of 1.5 m, were planted in the 52 m section of A1 channel which was located about 125 m upstream from the downstream weir and 400 m downstream from the upstream weir. The patches were planted in two configurations: the four most upstream patches in a dense configuration of 22 trees per m 2 and the three downstream patches in a sparse configuration of 7.3 trees per m 2 . The willow saplings were allowed to grow at the site for 10 months before the experiments started in August 2015. The average height of rooted willows was 0.4 m when they were planted, with an initial trunk diameter of 1 cm. In addition, indigenous small scale rough herbaceous vegetation and grasses were allowed to grow on the bank-side slopes of the channels. The bed itself was relatively unvegetated and mobile. Herbaceous vegetation had also spontaneously developed in between the willows in the vegetation patches. Although the flume was mostly dry during this in which the vegetation grew, although periods of prolonged flow occurred when the experiments took place in some (other) part of the flume. This regime of mostly dry and occasional flooding is typical for many Korean streams, which are ephemeral in nature. After the flow experiment, vegetation samples were taken to measure the height, diameter and morphology of the willows. The bed level in the experimental area was measured using RIEGL LMS-Z390i terrestrial laser scanner.
The flow measurements were taken after uniform, steady flow was achieved. The water depth in the channel downstream was fixed to approximately 1.1 m. Detailed three-dimensional velocity measurements with an Sontek 16 MHz ADV (Acoustic Doppler Velocimeter) were performed in the middle and 4 m upstream of the first dense patch (at D3, Fig. 2 ). The ADV devices were positioned pointing downward from an overhanging bridge. To get a signal within the vegetation patch, the devices were slightly moved to the left or right when needed. During trial convergence testing a required measurement time of 100 s was found. Both the discharge and the depth-averaged flow velocity profile are derived from the ADV measurements. We gathered data from three cases: (i) flow measurements at the upstream cross-section D0, (ii) flow measurements in the centre of the first willow patch, D3, (iii) and finally flow measurements at D3, but with the undergrowth removed ( Fig. 2 ). We will refer to these cases as D0, D3a(willows and undergrowth) and D3b (only willows) respectively. For case D0 and D3a ADV measurements were carried out every 30 cm, for a total of 23 locations in Y direction. For D3b a lower resolution was used of 60 cm, for a total of 12 locations in Y direction. In the vertical direction, a measurement was carried out every 5 cm. For D3a a vertical resolution of 2.5 cm was used in the vegetation patch.

Numerical experiment setup
A digital twin of the experimental flume was constructed with the open source Delft3D 4 modelling system ( Lesser et al., 2004 ), 1 using a two-dimensional, regular, numerical computational grid of 15 m × 171 m. The size of individual grid cells was 0.5 m (in flow direction) by 0.375 m with bed levels defined in the centre of each cell. A constant discharge, determined from the flume ADV measurements, was imposed at the upstream boundary, while the downstream boundary condition was defined as a constant water level to ensure a downstream water depth of approximately 1.1 m. The model was initialised with a constant water level. Each simulation ran for 30 min of model time with a time step of 0.3 s. The run time was chosen such that the entire model achieved equilibrium condition by the last time step. The scale of the model, given the relatively short run time, small water depths, fine grid and small time step -required changing some default numerical and physical parameters. We used a smoothing time of 2 min and a threshold depth for flooding/drying of 1 cm. The horizontal eddy viscosity was held constant throughout the model at 10 −2 m 2 s −1 . To be able to capture the sharp velocity gradients around the patches, low viscosity values were used ( Vionnet et al., 2004;Verschoren et al., 2016 ). In the experimental area, the geometric data was based on the bed level measurements of the experimental setup. For regions not covered by the bed level measurements, a synthetic cross-section was defined from the design specifications of the channel.
Bed friction was defined using four distinct roughness classes, two for the willow patches (dense and sparse), one for the vegetated channel slopes and one for the mobile sand channel bed. Each individual class is assigned its own roughness formula and friction parameters. Every grid cell boundary was assigned one of these four classes. All friction formulas are expressed in terms of the Chézy friction coefficient C [ m 1∕2 s −1 ]. For the channel bed we use the Manning friction formula: . The Keulegan equation (also known as the Colebrook-White equation) was found to be best suited to reproduce the depth-averaged velocity on the channel slopes: with Nikuradse roughness height k s [m] and parameters 1 = 18 m 1/2 s −1 , 2 = 12 . The values of both n b and k s were determined with GLUE.
There are various models available for resolving vegetation friction. Here, we used the two-layer approach of Baptist et al. (2007) , which performs favourably against laboratory experiments ( Vargas-Luna et al., 2015 ), and is generalised as follows: where C b is the Chézy friction coefficient of the bed without vegetation, [ − ] the vegetation parameter, g [ ms −2 ] the gravitational acceleration, the von Karman constant, h d [m] the deflected vegetation canopy height and an indicator such that = 0 for emergent conditions ( h < h d ) and = 1 for submerged conditions ( h > h d ).
The dimensionless vegetation term models the contribution of vegetation to the total friction of the vegetated part of the cross-section.
In the original formulation of Baptist et al. (2007) , it is computed as A fixed value for h d was used, which was estimated from measurements. Furthermore, we assumed that the friction in the patch was dominated by the vegetation and that the bed term had a neglible impact. Therefore we chose a constant value of C b at 60 m 1/2 s −1 , which led to negligible addition to the total roughness. The vegetation parameter was considered unknown, and estimated with GLUE.

GLUE method
The set of unknown parameters in the numerical experiment was given by = { , , } , i.e. the friction coefficients of the channel bed and slope and the vegetation parameter. The problem considered here was how to choose the values for parameters in such a way, that the measured flow velocities were reproduced by the flow model. While an optimal set of values may be found using one of various optimisation strategies available in literature, Beven (2006) argued that multiple parameter sets may well be found that all produce acceptable results. This creates uncertainty regarding the parameter values thus obtained through inverse modelling. The added benefit of GLUE above regular calibration procedures is that this uncertainty is formally taken into account. To achieve this, GLUE uses Bayesian inference, in which (usually very limited) prior knowledge about the parameter values is updated given the probability that those parameter values resulted in model output that compared favourably with measurements. This is technically achieved through Bayes' theorem: where is the error between modelled and measured flow velocities and u measured flow velocities. The prior distribution p ( ) that encodes our prior assumptions of probable parameter values for parameter set . The posterior distribution p ( | u ) is our goal, because it expresses the probability of parameter values in after having seen the measured flow velocities u . The likelihood  ( | , ) is a function of the model ( , ) and observations ( u ), and represents the probability of u given and .
To determine the functional form of the likelihood we needed to assume a statistical relationship between observed and modelled flow velocities. Here, we related the measured flow velocity vector u c at crosssection c to the simulated flow velocity ̂ , ( ) as follows: where c is a normal and independently distributed error term with variance 2 , i.e. ∼  (0 , 2 ) . The assumption of independent errors and zero bias expresses our assumption that the hydraulic model should be able to reproduce observe flow velocity profiles, such that any remaining discrepancy is adequately described by a normal distribution. As is commonly done in GLUE applications, we adopted a uniform prior distribution for each parameter. The upper and lower limits of each distribution were determined by exploratory computations with the model. For our adopted error model (5) , the corresponding likelihood function is: The unknowns in the inverse problem are both the model parameters and the variance of the residual errors 2 . To estimate 2 we used the maximum likelihood estimate 2 , = −1 ∑ =1 ( − ) 2 , with n the number of observations and the model results for the best performing parameter set ( Stedinger et al., 2008 ). Since the likelihood function returns the probability of the observed flow velocities given the model results evaluated with a given parameter set , model results that deviate significantly from observations will return likelihoods approaching zero. The likelihood function used here was formally derived from the adopted statistical model (5) . We note that GLUE allows for great flexibility in choosing from a variety of informal likelihood functions as well, which need an additional behavioural threshold to differentiate between behavioural and non-behavioural parameter sets. Since (6) tends to zero for very improbable parameter sets, such a behavioural threshold was not needed here.
The posterior parameter distributions, were obtained through Monte Carlo simulation using the Sobol' low discrepancy sequence ( Sobol', 1967 ). We sampled from the prior distributions to create a large number of possible parameter sets. Here, we used a sample size of 5000 model evaluations, with each evaluation using a different, randomly sampled set of parameter values. The same ensemble was used for cases D3a and D3b, using the value for h d derived from case D3a. Afterwards, the ensemble values of for D3b were corrected to account for a different height of the deflected vegetation canopy. The probability distribution for the parameter vector was obtained by application of Bayes' theorem. This procedure gave us two valuable sources of information on uncertainty related to model and observation. First, the estimate of 2 , or the 'predictive uncertainty'. This is the residual variance between model and measurement, which cannot be explained by choosing different parameter settings. The second is p ( | u ), i.e. the posterior parameter distribution or 'model uncertainty'. This expresses the uncertainty in the parameter values, given the measurements.

Flow measurements
Based on the ADV measurements, the flow fields covering the vertical (Y,Z) plane were constructed. To compute discharge from velocity measurements, it is usually necessary to assume a velocity profile (e.g. logarithmic, see Boiten, 2000 ). However, due to the irregular flow profiles expected in the flow through vegetation and the density of velocity measurements, we instead linearly interpolated points in between the ADV support points to a regular, cross-section covering grid, while assuming zero flow velocity at the bed ( Fig. 4 b, d, f). The total discharge was then computed by summing the product of the grid cell area. We found discharges of 2.91 m 3 s −1 (D0), 2.66 m 3 s −1 (D3a) and 2.54 m 3 s −1 (D3b). The differences are attributed to uncertainty in flow velocity measurements and interpolation inaccuracies.
The measurements at the unvegetated cross-section D0 show that the depth-averaged flow is low on the slopes and increases toward the centre of the channel, reaching about 0.5 ms −1 ( Fig. 4 a). At the right hand side of the channel ( y = 6 m) a local increase in flow velocity is measured. Given that the high velocities were consistently measured by different ADV devices, we assume they are not due to instrument problems. Therefore, we did not reject these measurements, but left it to the numerical model to explain whether these measurements are expected. We used the discharge at D3a as the upstream boundary condition in the numerical model. Two series of measurements were carried out for the cross-section in the willow patch (D3); one for the natural situation including organic debris and undergrowth (D3a) and the other for which the patch was cleaned (D3b).
Analysis of measured flow velocities for D3a revealed that one ADV measurement device returned near-zero results for all flow depths, suggesting malfunction. Measurements from this device (at = 6 . 4 m), were discarded from the results. Similar to D0, a region of higher flow velocities was observed at the right hand side of the patch near the water surface. In the patch itself, both the flow field ( Fig. 4 d) and the depth-averaged results ( Fig. 4 c) show slower velocities. However, the flow velocities near the water surface are similar to the unvegetated part of the channel, suggesting submerged flow. This is clear from the vertical velocity profile in the patch ( Fig. 3 ), which shows low flow ve-  locities ( < 0.3 ms −1 ) and large variation up to = −3 . 4 m, above which the velocities rapidly increase to more than 0.6 ms −1 and variation is significantly reduced. From Fig. 3 , we assumed that the vegetation was deflected such that the canopy height was 0.8 m. Case D3b -after removal of all organic debris and undergrowthbroadly shows ( Fig. 4 e, f) similar patterns to the D3a results, with some key differences in the flow through the patch. In the vertical profile ( Fig. 3 ) velocities near the bed ( < −3 . 6 m) are higher compared to D3a, which is attributed to removal of undergrowth. For higher water depths ( > −3 . 6 m) velocities show an initial decrease, while the variation over the patch increases. This is attributed to flow through the branches and leaves, which add comparatively more resistance. In contrast to D3a, clear submergence is not observed. Therefore we assume that in the case of D3b, the vegetation was just-submerged, meaning that the deflected canopy height was equal to the water depth.

Vegetation measurements
After the initial flow experiment (cases D0 and D3a), several samples were taken to determine the height and diameter of the willows ( Table 1 ). The stems were measured up until the upper knot, from which point several branches sprouted. We observed that organic debris which had attached itself around the willow stems increased the effective diameter of the plants by 50% on average. However, the distribution was not uniform -some plants experienced a significantly larger increase in diameter while others had little attached debris. This was reflected in the increased variance of the observations. The undergrowth had an average height of 29 cm and grew in between the willow stems ( Fig. 5 ) and consisted of various ways of herbaceous vegetation. It should be expected that both the undergrowth and the diameter increase due to organic debris increase the effective friction of the vegetation patch. The expected vegetation parameter is estimated based on the Baptist vegetation model ( = ℎ ) and the vegetation measurements from Table 1 . To compute , we assumed the branch height to be equal to the willow height minus the average stem height. Furthermore, we assumed rigid bending at the bed (following Verschoren et al., 2016 ), such that the deflected stem and branch heights can be computed from the deflected willow height (see Section 3.1 ) using standard trigonometry. We assumed a case drag coefficient of = 1 , which is a common assumption for submerged flow ( Wunder et al., 2011 ). Here, we did not account for the effect of foliage. The estimated vegetation parameter for case D3a, i.e. including debris, undergrowth and a deflected willow height of 0.8 m, is approximately 3 = 1 . 33 ± 1 . 044 . The large standard deviation in the expected parameter is mainly due to the variance in the diameter of the herbaceous vegetation in the undergrowth. For case D3b, i.e. no debris, no undergrowth and a deflected height equal to the willow height, is 3 = 0 . 76 ± 0 . 20 .

GLUE results
An important step within GLUE is to determine the ranges of the uniform prior distributions. These ranges should be generous, since it is important that these prior distributions cover the range where the (a priori unknown) posterior distributions will be. Assuming a trapezoidal channel, the lumped Manning coefficient for the given discharge and water depth is approximately 0.11 sm −1∕3 . Based on this relatively high friction factor and exploratory computations with the hydrodynamic model, we chose suitably large ranges for the prior distributions and a log-uniform distribution for ( Table 2 ).
For all three parameters ( n b , k s , ) we then computed the joint posterior distributions using GLUE. The posterior distributions give the likelihood that a certain parameter value leads to good model results. The width of those distributions is a measure of the parameter uncertainty regarding the possible 'true' value of the parameters given the model and the measurements. The marginal posterior distributions of all parameter values are shown in Fig. 6 .
The distribution of for D0 is very wide and does not show a clear peak ( Fig. 6 a). This indicates that the is an insensitive parameter for cross-section D0. Therefore, cannot be identified from flow measurements at that cross-section. This was not unexpected, as D0 is located 4 m upstream from the first patch and is therefore not directly influenced by the patch.
The other two distributions in Fig. 6 a show the values for that lead to good model results for D3a and D3b. D3a, which had significant undergrowth and debris, shows significantly larger inferred values for compared to D3b. The median value of 10 1.6 ( ~39.8) is an order of magnitude higher than would be expected based plant on parameters using the Baptist model. For the cleaned willows, case D3b, the values are much lower, with a median value of = 1 . 56 . This is more than twice as high as was expected based on the plant parameters. The uncertainty of the estimated parameter values for is significant, as can be observed by the width of the posterior distributions in Fig. 6 a and the standard deviation in Table 3 . For case D3a especially, the values for vary between 20 and 100, with some of the best simulations found near the upper boundary. It is important to note that at such high values, the sensitivity of to the roughness coefficient C is very much reduced due

Table 3
The number of observations n used for parameter estimation, the median and standard deviation of the model parameters in and the maximum likelihood estimation of . to the inverse square root in (3) , which may partly explain why such high uncertainty is found for this case. Therefore, values higher than the upper boundary are not expected to meaningfully improve results. Estimation of Manning coefficient of the channel bed n b and Nikuradse roughness height of the slopes k s resulted in well-defined posterior distributions for all cases ( Fig. 6 b,c), with standard deviations around 0.01 sm −1∕3 ( Table 3 ).
Interestingly, the distributions are not the same for the three cases. To measure the roughness of vegetation from water slope measurements, the roughness of the bed is commonly assumed to be independent from the roughness of the vegetation (e.g., Verschoren et al., 2016 ). However, the results shown in Fig. 6 suggest that the parameter values for both n b and k s are affected by the vegetation patch. For both n b and k s , the presence of a patch (D3a, D3b) results in higher parameter values. Therefore, the assumption of independence between the parameter of the vegetation patch ( ) and those of the channel bed and slope ( n b and k s ) would not have been valid in our case.
Modelled flow velocities, given the posterior probability distributions of the parameters, are compared with the depth-averaged ADV measurements ( Fig. 7 ). The depicted uncertainty bands show the range of the model uncertainty. Variation within the model uncertainty can be explained from uncertainty in the parameter values. Measurements that fall outside of these bands are explained through the residual error term in Eq. (5) . From Fig. 7 a we observe that the region of higher flow velocity in cross-section D0 cannot be explained by the numerical model. This indicates that this must be caused, either by an unmodelled process or feature, or by measurement error. However, in general the velocity profiles for all cases are well explained by the numerical model, which is reflected in the small standard deviation of the residual error ( , MLE ) for all three cases ( Table 3 ).
Finally, the uncertainty bands for cases D0 and D3a are comparatively small compared to the other cases, both in Figs. 6 and 7 . This is due to a higher density of ADV measurements. In general, a larger number of measurements helps to decrease model uncertainty and increase the identifiability of the individual parameters.

Identifiability of friction parameters
In this case study we identified three unknown friction parameters, related to the channel bed, vegetated channel slopes and the willow patch. Using only water level measurements, it is not possible to uniquely determine the values of these parameters. Our findings show that detailed measurements of the transverse depth-averaged velocity profile, in combination with a probabilistic inverse modelling approach (here, GLUE), allows us to estimate the parameter values and quantify the uncertainty of those estimations. The inverse modelling approach results in two different types of uncertainty. The first is model uncertainty, which is the uncertainty of ̂ , ( ) in (5) . This uncertainty can be decreased by increasing the number of observations (i.e. n in Table 3 ), or alternatively if a higher uncertainty is acceptable the number of observations may be decreased. In general, "data provides information and more data provides more information " ( Stedinger et al., 2008 ). The second type of uncertainty is predictive uncertainty, given by c in (5) . This uncertainty cannot be decreased without changing models, only more precisely estimated.
In literature it is often assumed that the total roughness is constituted of a linear sum of independent constituent terms (e.g. bed roughness and vegetation roughness). The total roughness can be estimated without detailed flow measurements. Therefore, if the bed roughness is known (e.g. from repeating the experiment without vegetation) the vegetation roughness can be determined. However, our results show that the assumption of independence between the terms would not have been valid for our study case. For both vegetated cases, the estimated friction of the bed and the slope was higher compared to the unvegetated case. A potential explanation for this could be that turbulent and shear stresses not sufficiently modelled by the parameterisation of eddy viscosity are compensated by a higher bed roughness.
The estimated parameter values for were found to be generally higher than was estimated based on the "rigid cylinder " estimator. This is especially true for case D3a, where the presence of undergrowth contributes to the overall failure of the estimator. In the cleaned case, D3b, the values are within the same order of magnitude, although still a factor two higher compared to the rigid cylinder estimation. A partial explanation could be found in the presence of foliage, which in other studies is reported to account for 60% of total drag ( Bal et al., 2011 ), as well as unmodelled stresses in the horizontal exhange layer between vegetated and unvegetated flow ( Truong et al., 2019 ).

Uncertainty of vegetation model parameters under field conditions
In this study we used the well-known two-layer model of Baptist et al. (2007) for submerged flow. While this model consistently Fig. 7. The modelled velocities (colored bands) and depth-averaged ADV measurements (markers) for the three cases. The colored uncertainty bands convey the result of the model uncertainty, while the dashed lines show the total predictive uncertainty. compares favourably to experimental data, application to field conditions is met with several challenges. An important limitation of the Baptist model is that reconfiguration is not explicitly taken into account. Reconfiguration can be taken into account in the two-layer approach by modelling two different effects. The first effect is streamlining of stems and foliage to decrease the total drag of the vegetation. This directly effects the vegetation parameter . The model proposed by Järvelä (2004) can be used to account for this. In this model, the vegetation parameter is defined by = The second effect of reconfiguration is deflection of the stems, which effectively decreases the height of the vegetation canopy height h d . This directly modifies the second term in the Baptist model which accounts for submerged flow and, depending on the chosen model, . For submerged flow, the estimation of vegetation friction is sensitive to the decrease of the vegetation canopy. This is especially so for low submergence ratios ( < 0.75) and high values of ( > 0.2), for which the contribution of vegetation friction (( 2 ) −1∕2 ) and contribution of the velocity profile above the vegetation ( √ ln ( ℎ ℎ )) to the total friction are of the same order of magnitude. Under these conditions, assuming rigid vegetation may lead to large errors in the estimation of vegetative friction.
The primary effect of undergrowth is an increased blockage ratio, which can be modelled by a higher value for . The inherent weakness of the rigid cylinder estimators for (used by Klopstra et al., 1997;Baptist et al., 2007;Huthoff et al., 2007;Yang and Choi, 2010 ) is that inclusion of undergrowth is contrived, while neglecting undergrowth may lead to significant underestimation of the effective roughness. Frontal blockage area estimators for (e.g. Västilä and Järvelä, 2017 ) do not suffer from this limitation, as undergrowth would be one part of the total blockage area. Comparing the two vegetated cases D3a and D3b, we observed lower flow velocities and greater deflection of the stems. This may be seen as a secondary effect of undergrowth, i.e. an increase in to-tal drag of the patch and therefore greater reduction in canopy height. This suggests that modelling this type of reconfiguration does not only depend on plant-specific parameters, but on the configuration of the entire patch, including secondary vegetation. Relatively simple predictive models based on single-plant behaviour, or even lumped models with a velocity dependent relationship with plant-specific parameters ( Verschoren et al., 2016 ), are therefore likely to include uncertainty regarding non-plant specific drivers.
Future effort may be directed to test an extended Baptist (or similar two-layer) model which would include the effects of streamlining and deflection under natural conditions. Ideally, these flume tests should cover a series of different discharges to test under varying water depths and flow velocities. To promote identifiability, the number of parameters of the extended model should be limited. Another argument for vegetative friction models with a relatively small number of parameters, apart from the issue of identifiability, is found in the usually datalimited problems in practice. Vegetation models with few parameters for which the uncertainty can be quantified are perhaps better suited for larger scale field applications than complex models which require intimate knowledge of vegetation configuration.

Conclusions
The objective of this study was to investigate flow resistance under natural vegetated conditions, in which secondary vegetation grows under and between the dominant species. Results show that the presence of undergrowth sharply increases the vegetation parameter , as well as increase the deflection of the willows. Overall, the vegetation parameter value was found to be higher than expected based on a priori estimations of a rigid cylinder estimator. This is attributed to aspects of real vegetation that deviate from the rigid-stick based approach, namely foliage, undergrowth and reconfiguration. This paper shows that velocity measurements in natural channels can be used effectively to estimate the parameter values of multiple sources of uncertainty under natural conditions using probabilistic inverse modelling. Results show well-defined posterior parameter distributions and model results compare favourably to measurements. The parameter values found in the non-vegetated cross-section were found to differ from the vegetated cross-section, which suggests that for the given model, the presence of a vegetation patch requires higher friction values in the surrounding non-vegetated part as well. An explanation for this has not yet been found.
Future challenges for practical application of vegetation models in large scale, 2D applications may require relatively simple models in which the presence of undergrowth and the effects of reconfiguration are taken into account. Inverse modelling of full scale experiments has been shown to be an appropriate tool for providing the experimental basis for field validation of these models.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.