# A topological restricted maximum likelihood (TopREML) approach to regionalize trended runoff signatures in stream networks

## Summary (6 min read)

### 1 Introduction

- Regionalizing runoff and streamflow for the purposes of making predictions in ungauged basins (PUB) continues to be one of the major contemporary challenges in hydrology.
- At global, regional and local scales only a small fraction of catchments are monitored for streamflow (Blöschl et al., 2013), and this fraction is at risk of decreasing given the ongoing challenge of maintaining existing gauging stations (Stokstad, 1999).
- Reliable information about local streamflows is essential for the management of water resources, especially in the context of changing climate, ecosystem and demography; flow prediction uncertainties are bound to propagate and lead to significantly suboptimal design and management decisions (e.g., Sivapalan et al., 2014).
- Techniques are needed to maximize the use of available data in datascarce regions to accurately predict streamflow, while providing a reliable estimate of the related modeling uncertainty.

### 1.1 Linear models

- There are a number of approaches to predicting runoff in ungauged catchments, including process-based modeling (e.g., Müller et al., 2014), graphical methods based on the construction of isolines (e.g., Bishop and Church, 1992), and statistical approaches.
- Such linear models are well understood and widely implemented, not only for PUB (see review in Blöschl et al., 2013, p.83) but also across a wide variety of fields in the physical and social sciences (e.g., Myers, 1990).
- Spatial correlation is generally problematic for linear model predictions, including the multiple regression approaches commonly used for regionalization.
- If the residuals are independent and identically distributed (iid), the best linear unbiased predictions (BLUP) of both y and its uncertainty (i.e., Var(y)) can readily be obtained using ordinary least squares (OLS) regression.
- Unfortunately, residuals are rarely iid in hydrological applications due to the spatial organization of hydrological processes around the topology of river channel networks.

### 1.2 Spatial correlation models

- There are several techniques available to address spatially correlated data.
- Huang and Yang, 1998), the theoretical justification for this approach is less robust than for point-scale processes.
- Runoff is organized around a topological network of stream channels, and the covariance structure implied for prediction should reflect the higher correlation between streamflow at watersheds that are “flow connected” (i.e., share one or more subcatchments), compared to unconnected but spatially proximate catchments.
- The bias increases quadratically with distance, meaning that estimates of the long-range variance (the sill) are strongly impacted by the presence of the trend, leading to an underestimation of the prediction uncertainty.
- The second type of approach embraces this topological structure.

### 1.3 The topological restricted maximum likelihood approach

- Rather than using a kriging estimator, the authors adopt a REML framework (Gilmour et al., 2004; Patterson and Thompson, 1971; Lark et al., 2006) to estimate variance parameters.
- This use of a REML framework to estimate a linear mixed effect model on a topological support is termed topological restricted maximum likelihood .
- For the local effects to form a suit- able basis for spatial interpolation, variations associated with temporal correlation (e.g., travel time effects) need to be removed, also known as – Treatment of time.
- First, TopREML is only suitable for the regionalization of time-averaged and statistically stationary runoff properties (i.e., runoff signatures).
- The underlying assumption is that runoff signatures of local flow generation regions that are close to each other (in Euclidian space) are more likely to be identical.

### 1.4 Paper outline

- The authors first derive the TopREML estimator and its variance for mass conserving (i.e., linearly aggregated) variables, with extensions to some non-conservative variables (Sect. 2).
- The authors then apply the approach in two case studies to evaluate its ability to predict mean runoff and runoff frequency by comparison to other available interpolation techniques: Sects. 3.1 and 4.1 present leave-one-out cross-validations in Nepal www.hydrol-earth-syst-sci.net/19/2925/2015/.
- In both cases, TopREML performed similarly to the best alternative geostatistical method.
- The authors then use numerical simulations to illustrate the effect of the two distinguishing features of TopREML: its ability to properly predict runoff using highly nested networks of stream gauges and its ability to properly estimate the prediction variance when accounting for observable features (Sects. 3.2 and 4.2).

### 2.2 Covariance structure of mass conserving variables

- In the linear mixed model framework, the propagation of hydrological variables through the flow network introduces Hydrol.
- Firstly, linearly propagated variables, such as annual specific runoff, are discussed.
- To account for the network structure, the catchment at any location along a stream can be subdivided into the IDA that are monitored for the first time by an upstream gauge.
- The authors assume that this function is well approximated by an exponential function ρ(ckm,φ)= exp(−c/φ).
- The nugget consists of the variance of processes that are spatially correlated over scales smaller than the sub-catchments and of measurement errors at the gauges.

### 2.4 E-BLUP and prediction variance at ungauged catchments

- Knowing φ̂, gn can be readily obtained from the relative position of site n and the gauges in the river network.
- The variance of the TopREML prediction error can be ex- pressed as Var(ỹn− yn)=.
- If is an error that is truly iid and does not affect the true value of yn (e.g., measurement errors), then Eq. (14) corresponds to the mean square error of the TopREML prediction of yn.
- If, by contrast, represents random variations of the true value of yn that are correlated over short distances (and so do not appear correlated in their data), then should be included in Eq. (11) and the prediction variance becomes Var(ỹn− yn)+ = Var(ỹn− yn)−+ σ 2, (15) because and u are independent.
- In reality is likely composed of both spatially correlated and iid error components and the true variance will be somewhere between these two bounds (Lark et al., 2006).

### 2.5 Application to non-conservative variables

- Unlike mean specific runoff, numerous streamflow signatures (e.g., runoff frequency or descriptors of the recession behavior) are non-conservative and cannot be expressed as linear combinations of their values in upstream subcatchments.
- In such conditions the derivations in Sect. 2.2 cannot be applied and the correlation structure in Eq. (7) will lead to biased REML predictions.
- Therefore, runoff frequency does not scale linearly through the river network.
- Therefore, although recession constants themselves do not propagate linearly, their value in ungauged basins can be estimated by taking the inverse of TopREML predictions of average response times.

### 2.6 Implementation

- To run the script, two vector data sets (e.g., ESRI Shapefiles) are needed as inputs – one containing the catchments where runoff is available and another containing the basins where predictions are to be made.
- Catchment polygons and explanatory and predicted variables must be provided as attributes of the vector polygons.
- The way in which the catchment polygons are nested provides the topology of the stream network.
- TopREML uses the Broyden–Fletcher–Goldfarb– Shanno (BFGS) algorithm (Wright and Nocedal, 1999) to maximize the restricted log likelihood, though stochastic algorithms are required if a non-differentiable (e.g., spherical) covariance function is selected.
- The authors found that initial values of [σ 20 ,φ0,ξ0] = [σ 2 LM,Eckm,1] worked well in their case studies, with σ 2LM the variance of the OLS residuals of the linear model and Eckm the average distance between IDA centroids.

### 3.1 Case studies

- Observed streamflow data are used to evaluate the ability of TopREML to predict streamflow signatures in ungauged basins.
- The assessment is based on leave-one-out cross-validations, where the tested model is applied to predict runoff at one basin based on observations from all the other basins.
- The Austrian data set was directly taken from the rtop package (Skøien et al., 2014), where mean summer runoff observations are provided to demonstrate Top-kriging.
- The Austrian data set did not contain additional observable features and previous studies have found spatial proximity to be a significantly better predictor of runoff than catchment attributes in Austria (Merz and Blöschl, 2005).

### 3.2.1 Network effects

- Conventional geostatistical methods predict runoff by weighing observations from surrounding basins based on their geographic distance.
- TopREML also incorporates the topology of the stream network by including or excluding basins based on their flow-connectedness.
- This adds topological information to the determination of the covariance structure of runoff, at the expense of discarding information that could be derived from correlations between spatially proximate regions that are not connected to the gauge of interest by a flow path.
- Assessing the net benefits of accounting for network effects requires being able to control the topology of the network, and thus requires numerical simulations.
- A non-topological geostatistical method like universal kriging would include all basins within and exclude all basins beyond the spatial autocorrelation range.

### 3.2.2 Variance estimation and observable features

- A key advantage of the REML framework is its ability to avoid the downward bias in the covariance function that affects kriging-based methods (including Top-kriging) when external trend coefficients are simultaneously estimated.
- Again, empirical cross-validation analysis does not allow an assessment of this bias, because the observation data sets used contained only one observation per location.
- We construct the observed prediction uncertainty by taking the standard deviation of the prediction errors across all 1000 Monte Carlo runs and compare it to the square root of the median predicted variance.the authors.the authors.
- The authors compare TopREML and Top-kriging based on their ability to model prediction variance.
- Conversely, including a trend in the model will cause the remaining error to mostly consists of (spatially uncorrelated) residuals, so in this case Eq. (14) is used.

### 4.1 Case studies

- Basin-level predictions of the considered signatures are presented in Fig. 5 for the three cross-validation analyses described in Sect. 3.1. Figure 5 also provides box plots summarizing the distribution of the ensuing cross-validation errors.
- In the three analyses, the prediction errors related to TopREML were comparable to the best alternate method: a linear model for annual specific runoff and Topkriging for runoff frequency and summer runoff .
- Allowing for spatial correlation in the residuals (UK) decreased the median absolute error by 11 % compared to the linear model (LM) for runoff frequency in Nepal and 31 % for summer runoff in Austria.

### 4.2 Numerical simulation

- Results from the Monte Carlo analysis are presented in Fig. 6, showing the outcomes of the two numerical experiments described in Sect. 3.2. Figure 6a and b shows the effect of network complexity on the performance of TopREML relative to the baseline performance of universal kriging.
- The effect is expected to increase with Nouter and decrease with Ninner, reaching zero when 100 % of observed basins lie within the spatial correlation range and 0 % of the basins beyond the range are flow-connected.
- A linear regression of the relative performance of TopREML against Nouter and Ninner showed that both trends are significant and in the expected direction.
- In Fig. 6c, the Monte Carlo analysis showed that model uncertainty is well predicted by TopREML and strongly underestimated by Top-kriging, both with and without considering an external trend.

### 5.1 Performance of TopREML

- Cross-validation outcomes suggest that TopREML is an attractive operational tool for predicting streamflow in ungauged basins.
- Both results suggest that the benefit of accounting for network effects on correlations outweighs the cost of losing some information on the correlation between unconnected basins.
- This allows TopREML to accurately predict modeling uncertainties even for highly trended and auto-correlated runoff signatures, as visible in the Monte Carlo analysis presented on Fig. 6c.
- TopREML also has considerably lower computational requirements than Top-kriging, both in terms of input data and optimization complexity.
- Indeed, TopREML does not rely on a distributed point process but assumes homogenous IDAs.

### 5.2 Model selection

- The regionalization methods assessed in the cross-validation analysis range from simple linear regressions with strong independence assumptions, to complex geostatistical methods that allow for both spatial and topological correlations.
- In that situation, parsimony prescribes selecting the least complex of the best performing methods.
- A dense network of flow gauges is necessary for geostatistical methods to properly estimate the semivariogram and improve on predictions from linear regressions – the case studies suggest that the mean distance between the gauges must be on the order of half the spatial correlation range of the runoff signature.
- Including these observable features in the model reduces the correlation scale of the residuals, possibly crossing the threshold below which geostatistics are not the most parsimonious approach.
- In that case there is a tradeoff between relying on observable features or variance information to make a prediction, and parsimony and stationarity considerations come into play when selecting the regionalization model.

### 6 Conclusions

- The approach takes into account the spatially correlated nature of runoff and the nested character of streamflow networks.
- The method was successfully tested in cross-validation analyses on mass conserving (mean streamflow) and nonconservative (runoff frequency) runoff signatures in Nepal (sparsely gauged) and Austria (densely gauged), where it matched the performance of the best alternative method: Topkriging in Austria and linear regression in Nepal.
- TopREML outperformed Top-kriging in the prediction of uncertainty in Monte Carlo simulations and its performance is robust to the inclusion of observable features.
- This flexibility, along with its ability to provide a reliable estimate of the prediction uncertainty, offer considerable scope to use this computationally inexpensive method for practical PUB applications.

Did you find this useful? Give us your feedback

##### Citations

237 citations

24 citations

19 citations

19 citations

### Cites background from "A topological restricted maximum li..."

...Possibly, the most obvious stream variable is runoff (Müller and Thompson, 2015), however, since there is a correlation between flow-rate and the dilution capacity of streams, there aremanyother variables related towater-flow, such as the measurement of water physicochemical variables,…...

[...]

...New models of covariance valid for stream networks (Laaha et al., 2012; Müller and Thompson, 2015) v....

[...]

16 citations

##### References

5,555 citations

4,487 citations

3,855 citations

### "A topological restricted maximum li..." refers methods in this paper

...Rather than using a kriging estimator, we adopt a Restricted Maximum Likelihood (REML) framework (Gilmour et al., 2004; Patterson and Thompson, 1971; Lark et al., 2006) to estimate variance parameters....

[...]

2,500 citations

### "A topological restricted maximum li..." refers background or methods in this paper

...…from a second order stationary ran- dom process, and that the covariance between y0 k and y0 m will depend only on the relative position of sub-catchments m and k, given some specified correlation function ⇢(·) of the distance c km between the centroids of the two sub catchments (Cressie, 1993)....

[...]

...This reduces the bias on the semivariogram by allowing the variance to be estimated independently from the trend120 coefficients (Cressie, 1993; Lark et al., 2006)....

[...]

...By contrast, the expected downward bias in the kriging estimation of partial sills (Cressie, 1993) is clearly visible in the underestimation of prediction uncertainties by the Top-kriging method....

[...]

...Within PUB, kriging (Cressie, 1993) based geostatistical methods have been widely used (e.g, Huang and Yang, 1998; Gottschalk et al., 2006; Sauquet, 2006; Sauquet et al., 2000; Skøien et al., 2006)....

[...]

...Following Cressie (1993) (p. 68), the covariance between two aggregated random variables y0 k and y0 m is expressed as a function of the covariogram C P (·) of the underlying point-scale process: Cov (y0 k ,y0 m ) = 1 A2 Z Sk Z Sm C P (| x2 x1 |)dx1dx2 = 1Z 0 ⌫(D)C P (D)dD (A1)575 where S k and S m…...

[...]

2,455 citations