Inferred rainfall sequences generated by a novel method of inverting a continuous time transfer function show a smoothed profile when compared to the observed rainfall, however, streamflow generated using the inferred catchment rainfall is almost identical to observed streamflow (Rt2 > 97%). This paper compares the effective rainfall inferred by the regularised inversion process (termed inferred effective rainfall (IER)) proposed by the authors with effective rainfall derived from the observed catchment rainfall (termed observed effective rainfall (OER)) in both time and frequency domains in order to confirm that, by using the dominant catchment dynamics in the inversion process, the main characteristics of catchment rainfall are being captured by the IER estimates. Estimates of the resolution of the IER are found in the time domain by comparison with aggregated sequences of OER, and in the frequency domain by comparing the amplitude spectra of observed and IER. The temporal resolution of the rainfall estimates is affected by the slow time constant of the catchment, reflecting the presence of slow hydrological pathways, for example, aquifers, and by the rainfall regime, for example, dominance of convective or frontal rainfall. It is also affected by the goodness-of-fit of the original forward rainfall–streamflow model.

INTRODUCTION

Rainfall is the key driver of catchment processes and is usually the main input to rainfall–streamflow models. If the rainfall and/or streamflow data used to identify or calibrate a model are wrong or disinformative, the model will be wrong and cannot be used to predict the future with any certainty. Blöschl et al. (2013) state that if the dominant pathways, storage and time-scales of a catchment are well defined then a model should potentially reproduce the catchment dynamics under a range of conditions. It is often the case that hydrological variables, such as rainfall and streamflow, are measured at hourly or sub-hourly intervals then aggregated up to a coarser resolution before being used as input to rainfall–streamflow models resulting in the loss of information about the finer detail of the catchment processes (Littlewood & Croke 2008, 2013; Littlewood et al. 2010). Kretzschmar et al. (2014) have proposed a method for inferring catchment rainfall using sub-hourly streamflow data. The resulting rainfall record is smoothed to a coarser resolution than the original data but should still retain the most pertinent information.

This paper investigates the implications of the reduced resolution and the potential loss of information introduced by the regularisation process in both the time and frequency domains. Both temporal and spatial aggregation are incorporated in the transfer function model, however, only the temporal aspect is considered here. The effect of spatial rainfall distribution using sub-catchments will be the subject of a future publication.

The method developed and tested by Kretzschmar et al. (2014) – termed the RegDer method – inverts a continuous time transfer function (CT-TF) model using a regularised derivative technique to infer catchment effective rainfall from streamflow with the aim of improving estimates of catchment rainfall arguing that a model that is well-fitting and invertible is likely to be robust in terms of replicating the catchment system. In the context of this study, observed catchment rainfall (may be derived from one or more rain-gauges by any suitable method, e.g., Thiessen polygons) is converted to observed effective rainfall (OER) by a non-linear transform designed to render the relationship between the rainfall input and streamflow output (via a CT-TF) linear. The inversion process takes the catchment streamflow and, using a regularisation process, infers effective rainfall (IER), which is then converted to inferred catchment rainfall (ICR) by the reverse of the non-linear transform.

The effective rainfall (both OER and IER) may be termed scaled rainfall (related to Andrews et al. (2010)) as it is derived from the overall catchment rainfall.

The classical approach to inverse (as opposed to reverse) modelling involves the estimation of non-linearity (rainfall or baseflow separation) and the unit hydrograph (UH), which is an approximation to the impulse response of the catchment. Boorman (1989) and Chapman (1996) use sets of event hydrographs to estimate the catchment UH. Boorman (1989) superimposed event data before applying a separation technique and concluded that the data required may be more coarsely sampled than might be expected because one rain-gauge is unlikely to be representative of the whole catchment. Chapman (1996) used an iterative procedure to infer rainfall patterns for individual events before applying baseflow separation. The resultant UHs had higher peaks and shorter rise times and durations than those obtained by conventional methods. He viewed the effective rainfall as the output from a non-linear store. Duband et al. (1993) and Olivera & Maidment (1999) used deconvolution to identify mean catchment effective rainfall, which was redistributed using relative runoff coefficients while Young & Beven (1994) based a method for inferring effective rainfall patterns on the identification of a linear transfer function. In that study, a gain parameter varying with time accounted for the non-linearity in the relationship between rainfall and streamflow.

In recent years, a range of different approaches has been used to explore reverse modelling in hydrology, that is, estimating effective rainfall from streamflow. Notable publications include Croke (2006, 2010), Kirchner (2009), Andrews et al. (2010), Young & Sumisławska (2012), Brocca et al. (2013, 2014) and Kretzschmar et al. (2014). Kirchner's method links rainfall, evapotranspiration and streamflow through a sensitivity function making assumptions which allow rainfall to be inferred from the catchment streamflow. The method has been applied by Teuling et al. (2010) and Krier et al. (2012) to catchments in Switzerland and Luxembourg and has been found to work for catchments with simple storage–streamflow relationships and limited hysteresis. Brocca et al. (2013) employed a similar method based on the water balance equation but inferred the rainfall series from soil moisture. In a further study, Brocca et al. (2014) used satellite-derived soil moisture to infer global rainfall estimates. Rusjan & Mikoš (2015) applied Kirchner's simple dynamic system concept to a catchment in south-west Slovenia characterised most of the time by subsurface storage but showing a response that by-passed this storage after intense rainfall. They combined two separate sensitivity functions to enable the simulation of a range of contrasting hydrological conditions. Croke (2006) derived an event-based UH from streamflow alone but his approach was limited to ephemeral quick-flow-dominant catchments while Andrews et al. (2010) and Young & Sumisławska (2012) used a discrete model formulation inverted directly or via a feedback model (which could be adapted to CT formulation). Croke (2010) explored a similar approach to the one presented in Andrews et al. (2010) for several catchments. This is done in the context of slow flow, recharge and quickflow separation, relating the derived general model to existing ones (such as IHACRES). He also includes measures to constrain the rainfall estimate uncertainty. The flow components are estimated as individual discrete-time transfer functions separated using a relaxation procedure. The equivalent effective rainfall estimate is then obtained as a form of inverse discrete-time transfer function with the separated flow components as inputs. The approach proposed by Kretzschmar et al. (2014) combined a CT-TF model with regularized derivative estimates to infer the catchment rainfall from sub-hourly streamflow data, including comparisons to the direct inverse of a discrete transfer function model, similar to those used by Croke (2010) and Andrews et al. (2010).

Littlewood (2007) applied the IHACRES model (e.g., Jakeman et al. 1990) to the River Wye gauged at Cefn Brwyn showing that the values for the model parameters for that catchment changed substantially as the data time-step used for model calibration decreased. Littlewood & Croke (2008) extended this work to include a second catchment and found that as the time-step decreased the parameter values approached an asymptotic level (on a semi-log plot) concluding that, at small enough time-steps, parameters become independent of the sampling interval. They suggested further investigation using data-based mechanistic (DBM) modelling methods as described by Young & Romanowicz (2004) and Young & Garnier (2006) for estimating CT models from discrete input data. Such models generate parameter values independent of the input sampling rate – as long as the sampling rate is sufficiently high in comparison to the dominant dynamics of the system. Advantages of using the CT formulation include allowing a much larger range of system dynamics to be modelled, e.g., ‘stiff’ systems that have a wide range of time constants (TC), typical of many hydrological systems. The outputs from such a model can be sampled at any time-step, including non-integer, and the parameters have a direct physical interpretation (Young 2010).

Krajewski et al. (1991) compared the results from a semi-distributed model and a lumped model and concluded that catchment response is more sensitive to rainfall resolution in time than space, while a study by Holman-Dodds et al. (1999) demonstrated that models calibrated using a smoothed rainfall signal (due to coarse sampling) may result in underestimation of streamflow. Further calibration, required to compensate, leads to the loss of physical meaning of parameters. They also concluded that parameters estimated at one sampling interval were not transferable to other intervals; a conclusion echoed by Littlewood (2007) and Littlewood & Croke (2008).

Studies by Clark & Kavetski (2010) showed that in some cases, numerical errors due to the time-step are larger than model structural errors and can even balance them out to produce good results. The follow-up study by Kavetski & Clark (2010) looked at its impact on sensitivity analysis, parameter optimisation and Monte Carlo uncertainty analysis. They concluded that use of an inappropriate time-step can lead to erroneous and inconsistent estimates of model parameters and obscure the identification of hydrological processes and catchment behaviour. Littlewood & Croke (2013) found that a discrete model using daily data overestimated TC for the River Wye gauged at Cefn Brwyn when compared to those estimated from hourly data confirming that parameter values were dependent on the time-step. They discussed the loss of information due to the effect of time-step on time constants and suggested that plots of parameter values against time-step could be used as a model assessment tool. In a previous study, Littlewood & Croke (2008) compared the sensitivity of parameters for two catchments with respect of time-step and discussed the role of time-step dependency on the reduction of uncertainty. They also suggested CT-TF modelling using sub-hourly data to derive sampling rate independent parameter values. Littlewood et al. (2010) introduced the concept of the Nyquist–Shannon (N–S) sampling theorem, which defines the upper bound on the size of sampling interval required to identify the CT signal without aliasing, and consequentially its effect on the frequency of sampling required to specify a rainfall–streamflow model. Given a frequent enough sampling rate, the CT model is time independent and can be interpreted at any interval.

Further understanding may be gained by transforming rainfall and streamflow series from the time domain to the frequency domain and using spectral analysis. Several potential uses of spectral analysis in hydrology have been explored including modelling ungauged catchments, modelling karst systems and seasonal adjustment of hydrological data series. A maximum likelihood method for model calibration based on the spectral density function (SDF) has been suggested by Montanari & Toth (2007). The SDF can be inferred from sparse historic records in the absence of other suitable data making it a potentially useful tool for modelling ungauged catchments. They also suggest that spectral analysis may provide a means of choosing between different apparently behavioural models. Cuchi et al. (2014) used ‘black box’ modelling and frequency analysis to study the behaviour of a karst system (located at Fuenmajor, Huesca, Spain). They concluded that the method works well for a linear system and that Fuenmajor has a linear hydrological response to rainfall at all except high frequencies. They suggest that the non-linearity issues might be addressed using appropriate techniques such as wavelets or neural networks. Szolgayová et al. (2014) utilised wavelets to deseasonalise a hydrological time-series and suggested that the technique had potential for modelling series showing long-term dependency (interpreted as containing low frequency components).

The method introduced by Kretzschmar et al. (2014) showed that given that the rainfall–streamflow model captures the dynamics of the catchment system, the high frequency detail of the rainfall distribution is not necessary for the prediction of streamflow due to the damping (or low-pass filter) effect of the catchment response. The numerical properties of the regularisation as applied to the inversion process place a mathematical constraint of smoothness balanced against a loss of some temporal resolution in the inferred rainfall time-series. The regularisation and therefore smoothing level is controlled through the Noise Variance Ratio (NVR), optimised as part of the process and is only applied when necessary, i.e., when the analytically inverted catchment transfer function model is improper (has a numerator order higher than the denominator order).

APPLICATION CATCHMENTS

RegDer has been tested on two headwater catchments with widely differing rainfall and response characteristics – Baru in humid, tropical Borneo and Blind Beck, in humid temperate UK.

Baru – tropical catchment

The 0.44 km2 Baru catchment (Figure 1(a)) is situated in the headwaters of the Segama river located in Sabah on the northern tip of Borneo, East Malaysia (4° 580′ N 117° 490′ E). The climate is equatorial with a 26 year (1985–2010) mean rainfall of 2,849 mm (Walsh et al. 2011) showing no marked seasonality but tending to fall in short (<15 min) convective events showing high spatial variability and intensities much higher than those of temperate UK (Bidin & Chappell 2003, 2006). Due to the high spatial variability, a network of six automatic rain-gauges (13.6 gauges per km2) was used to derive the catchment-average rainfall using the Thiessen polygon method. Haplic alisols, typically 1.5 m in depth and with a high infiltration capacity (Chappell et al. 1998), are underlain by relatively impermeable mudstone bedrock resulting in the dominance of comparatively shallow sub-surface pathways in this basin (Chappell et al. 2006, 2007). As a result of the high rainfall intensity and shallow water pathways, the stream response is very flashy (i.e., rapid recession in the impulse response function). Vegetation cover is lowland, evergreen dipterocarp forest, which was subject to selective logging during 1988–1989 (Greer et al. 1993). The data used in the analysis are from February 1996 sampled at 5 min intervals (Figure 1(b)) and have been modelled previously by Chappell et al. (1999) and Walsh et al. (2011).

Figure 1

(a) The location of the 0.44 km2 tropical Baru catchment in Sabah (dark grey area in bottom left map – Sabah Foundation Forest management concession), Borneo and (b) the hydro- and hyetographs for the February 1996 sampled at 5 min intervals showing the flashy response of the catchment to the high intensity, spatially variable rainfall.

Figure 1

(a) The location of the 0.44 km2 tropical Baru catchment in Sabah (dark grey area in bottom left map – Sabah Foundation Forest management concession), Borneo and (b) the hydro- and hyetographs for the February 1996 sampled at 5 min intervals showing the flashy response of the catchment to the high intensity, spatially variable rainfall.

Blind Beck – temperate catchment

The Blind Beck catchment (Figure 2(a)) has an area of 8.8 km2 and lies in the headwaters of the Eden basin in northwest England, UK (54.51 °N 2.38 °W). The basin's response shows evidence of deep hydrological pathways due to the presence of deep limestone (62%) and sandstone (38%) aquifers resulting in a damped hydrograph response (Mayes et al. 2006; Ockenden & Chappell 2011; Ockenden et al. 2014). Winter rainfall in this basin is derived from frontal systems with typically lower intensities than the convective systems in the tropics (Reynard & Stewart 1993). Data from a single tipping bucket rain-gauge (i.e., 0.1 gauges per km2) located in the middle of the catchment were used in this study. The data used in the analysis cover the period from 26th December 2007 at 16:45 to 31st December 2007 at 21:45 sampled at 15 min intervals (Figure 2(b)) and was previously modelled by Ockenden & Chappell (2011) using an aggregated hourly time-step.

Figure 2

(a) The location of the 8.8 km2 temperate Blind Beck catchment in northwest England and (b) the hydro- and hyetographs for Blind Beck for the period from 26th December 2007 at 16:45 to 31st December 2007 at 21:45 sampled at 15 min intervals showing its response to less intense frontal rainfall and deeper hydrological pathways.

Figure 2

(a) The location of the 8.8 km2 temperate Blind Beck catchment in northwest England and (b) the hydro- and hyetographs for Blind Beck for the period from 26th December 2007 at 16:45 to 31st December 2007 at 21:45 sampled at 15 min intervals showing its response to less intense frontal rainfall and deeper hydrological pathways.

The choice of these two experimental catchments, therefore, allowed the initial evaluation of the estimation of catchment rainfall from streamflow for the end-member extremes of a basin with tropical convective rainfall and shallow flow pathways to a basin with temperate frontal rainfall (i.e., much lower intensity) and deep flow pathways (i.e., much greater basin damping or temporal integration).

MODEL FORMULATION AND PHYSICAL INTERPRETATION

This study investigated the limits of inferred catchment effective rainfall estimation from streamflow. CT-TF models identified from the observed data using DBM modelling approaches (Young & Beven 1994; Young & Garnier 2006), are inverted using the RegDer method (Kretzschmar et al. 2014) and used to transform catchment streamflow into estimates of catchment inferred rainfall.

DBM modelling makes no prior assumptions about the model structure (although it often uses structures based on transfer functions), which is suggested by the observed data, and must be capable of physical interpretation. As transfer functions are linear operators, a transform structured as a bilinear power-law (Equation (1)), also identified from the observed data, was applied to linearise the data before model fitting (Young & Beven 1994; Beven 2012, p. 91): 
formula
1
where P is the observed rainfall, Q the observed streamflow and α is a parameter, estimated from the data. Pe is the effective observed rainfall (ERER) and Q is used as a surrogate for catchment wetness. Both catchments used in this study proved to be predominantly linear in their behaviour so transformation (Equation (1)) was not used. In the initial study, a wide range of possible models was identified using algorithms from the Captain Toolbox for Matlab (Taylor et al. 2007). The models selected were a good fit to the data and were suitable for inversion. The Nash–Sutcliffe efficiency (NSE or Rt2) is commonly used to compare the performance of hydrological models. Often, several models can be identified that fit the data well (the equifinality concept of Beven (2006)). From these, models with few parameters to be estimated that inverted well were selected. In this study, a second-order linear model was found to fit both catchments. The output from the RegDer process is an IER series to which the inverse of the power law is then applied, if necessary, to construct an ICR sequence. The process is illustrated in Figure 3.
Figure 3

Model identification and inversion workflow where P is the observed catchment rainfall, Pe is the effective rainfall, Q is the observed streamflow, Peh is the IER and Ph the ICR. Non-linearity is represented by the bilinear power law (Beven 2012, p. 91). The CT-TF is given by G(s) where A(s) and B(s) are the denominator and numerator polynomials and the inversion process is represented by G−1(s) where A*(s) and B*(s) refer to the symbolic denominator and numerator polynomials of the regularised inverse transfer function as in Equation (4).

Figure 3

Model identification and inversion workflow where P is the observed catchment rainfall, Pe is the effective rainfall, Q is the observed streamflow, Peh is the IER and Ph the ICR. Non-linearity is represented by the bilinear power law (Beven 2012, p. 91). The CT-TF is given by G(s) where A(s) and B(s) are the denominator and numerator polynomials and the inversion process is represented by G−1(s) where A*(s) and B*(s) refer to the symbolic denominator and numerator polynomials of the regularised inverse transfer function as in Equation (4).

The transfer function model inversion process has been described in Kretzschmar et al. (2014). It involves transition from the transfer function catchment model: 
formula
2
to the direct inverse (in general non-realisable): 
formula
3
which is then implemented using regularised streamflow derivatives in the form of: 
formula
4
where is the Laplace transform of the optimised regularised estimate of the nth time derivative of Q:. The regularised derivative estimates replace the higher order derivatives in Equation (3), which otherwise make Equation (3) unrealisable (improper) – this is the core of the method in Kretzschmar et al. (2014). In the implementation, nth derivatives in Equation (4) are not estimated, but advantage is taken of the filtering with the denominator polynomial, whereby only (n-m)th order regularised derivative estimates of Q are required in combination with a proper transfer function.

The IER sequences generated by RegDer generally have a much smoother profile (illustrated in Figure 4) than the OER inputs, however, streamflow sequences generated with the ICR used as the model input are almost indistinguishable from the observed streamflow (Rt2 > 97%). This indicates that the catchment dynamics, as captured by the transfer function model, renders the differences between observed and inferred rainfall immaterial. The reason for this becomes clear when looking at the frequency domain analysis of the inversion process shown in this paper.

Figure 4

Observed effective and IER profiles generated using the RegDer inversion method for (a) Blind Beck and (b) Baru.

Figure 4

Observed effective and IER profiles generated using the RegDer inversion method for (a) Blind Beck and (b) Baru.

In order to investigate this, the IER is compared to aggregated effective observed rainfall sequences with increasing levels of aggregation until a good match is found (high value of Rt2 or R). Two methods of aggregation have been used: (1) averaging over a range of time-series and (2) moving average over varying time-scales. Two measures are used to assess the correspondence between the IR and the aggregated effective rain: (1) Rt2, the coefficient of determination and (2) R, the instantaneous Pearson correlation coefficient. They are given by: 
formula
5a
 
formula
5b
where ER indicates a value from the aggregated effective rainfall sequence with mean and IER is the corresponding value from the IER sequence with mean . Both Rt2 and R values tend towards a maximum value as aggregation increases. The aggregation level at which the maximum is reached is identified and taken as an estimate of the resolution of the inferred effective series. This value is then compared to the system fast time constant (TCq) and the N–S sampling limit.

Continuous model formulation

One of the advantages of using a CT model formulation is that the parameters have a direct physical interpretation independent of the model's sampling rate (Young 2010). The continuous time model formulation for a second-order model is given by: 
formula
6
where y is the measured streamflow at time t, is the transport delay and u is the effective rainfall at time t. If the denominator can be factorised and has real roots, Equation (6) can be rewritten as: 
formula
7
where TCq and TCs are the system time constants and are often significantly different – a ‘stiff’ system. Decomposing the model into a parallel form gives: 
formula
8
where gq and TCq are the steady state gain and time constant of the fast response component and gs and TCs are the steady state gain and time constant of the slow response component. The steady state gain of the system as a whole is given by: 
formula
9
so the fraction of the total streamflow along each pathway can be calculated from: 
formula
10

The fraction of streamflow attributed to the slow response component is sometimes termed the slow flow index (SFI) (Littlewood et al. 2010). The example shown here uses a second-order model but the general principle can be extended to higher order models. Details of the inversion and regularisation processes can be found in Kretzschmar et al. (2014).

Sampling frequency

The N–S frequency gives the upper limit on the size of the sampling interval, Δt, that will enable the system dynamics to be represented without distortion (aliasing – Bloomfield 1976, p. 21). Aliasing occurs when a system is measured at an insufficient sampling rate to adequately define the signal from the data.

The N–S theorem states that the longest sampling step for a signal with bandwidth (maximum frequency, where in cycles per time unit) to be represented is: 
formula
11
in order to completely define the system in absence of observation disturbance (Young 2010). If the sampling interval is small enough to uniquely define the system, the estimated CT model should be independent of the rate of sampling. Conversely, if the frequency of the inferred output is less than the N–S limit, then the system dynamics should be adequately captured. Other estimates of the sufficient sampling interval, designed to avoid proximity to the Nyquist limit, have been made by Ljung (1999) and Young (2010). In terms of system TCs, these limits are given by: 
formula
12a
 
formula
12b
 
formula
12c

Temporal aggregation of effective rainfall

Two methods for aggregating ER were used to estimate the time resolution of the IER. Rainfall is the total volume accumulated over the sampling interval so the ER was aggregated over progressively longer sampling periods of 2 to 24 times the base sampling period and averaged to form a new smoothed sequence that could be compared with the IER. For comparison, aggregation was also performed via a moving average process utilising the convolution method available in Matlab. Both methods may be affected by the aggregation starting point and edge effects. The aggregated ER sequences were compared to the IER using the coefficient of determination (Rt2) and the correlation (R). Rt2 and R tend towards a maximum value as aggregation increases. The aggregation time-step at which this value is established is used to estimate the resolution of the IER.

Spectral analysis

Periodograms of the amplitude spectra of the observed and modelled series were plotted to test whether the ER and IER have the same dynamics in the critical frequency range, despite the loss of time resolution (related to low pass filtering due to regularisation). A periodogram is the frequency domain representation of a signal; transforming the signal into the frequency domain may reveal information that is not visible in the time domain. A transfer function shown in its equivalent frequency domain form describes the mapping between the input and the output signals' spectra for the linear dynamic systems used here. Signals may be easily transformed between the time and frequency domains (Wickert 2013).

Periodograms are obtained using the Matlab implementation of the fast Fourier transform and smoothed using the integrated random walk (e.g. Young et al. 1999); the same regularisation approach as used in the calculation of the IER, implemented in the Captain Toolbox (Taylor et al. 2007). Periodograms of ER, IER and catchment streamflow are compared on a single plot showing how the spectral properties of the inversion process are used to obtain the IER. The streamflow spectrum is the result of mapping the rainfall spectrum by the catchment dynamics. To make a full inversion of that mapping would involve very strong amplification of high frequencies with all the negative consequences discussed by Kretzschmar et al. (2014). The most significant implications of full inversion include the introduction of high amplitude, high frequency noise artefacts into the rainfall estimates. The regularisation of estimated derivatives introduces the effect of low-pass filtering into the inversion process, avoiding the excessive high frequency noise. Regularisation does not introduce any lag into the process, unlike traditional low-pass filtering.

RESULTS AND DISCUSSION

Figure 4 illustrates the smoothed rainfall distribution of the IER sequence obtained using the RegDer method. Similar streamflow sequences are generated using either the observed rainfall or ICR sequences as model input (see Kretzschmar et al. 2014). The implication is that the catchment system dynamics are being captured despite the apparent difference in the rainfall distribution and that the detail of the rainfall series in time may not be important when modelling the dominant mode of streamflow dynamics.

In order to assess the degree of resolution lost by estimating rainfall using the RegDer method, the ER was aggregated using two methods (i.e., simple aggregation by resampling and a moving average) and the resulting sequences compared to the IER sequence in the time domain. Plots of progressively more aggregated sequences are shown in Figure 5. It can be seen that as aggregation increases, peaks become lower and more spread out and the sequence is effectively smoothed. The coefficient of determination (Rt2) and the correlation (R) between the aggregated sequence and the IER tends to a maximum then decreases as aggregation time increases – ultimately the variation in the sequence would be completely smoothed out. The point at which the maximum value is reached is taken as an estimate of the resolution of the IER. Plots of Rt2 or R values are shown in Figures 6 (aggregation by resampling) and 7 (moving average estimate). Time resolution estimates are shown in Table 1 and compared with the fast time constant (TCq) and the N–S sampling limit.

Table 1

The sampling frequency (hr), time constants (TCq – fast and TCs – slow), SFI (the percentage of the flow taking the slow pathway), the N–S safe sampling limit (hr) and the time resolution of the IER estimated by resampling and moving average methods. Also shown is the frequency domain estimate of the resolution – the cut-off point beyond which the signal carries very little information (illustrated in Figure 8) and can be considered unimportant. The estimated time resolution of the IER sequences is less than the dominant (fast) mode of the catchments and considerably less than the ‘safe’ N–S limit

            Time resolution estimates
 
  
Catchment Sampling frequency (hr) TCq (hr) TCs (hr) SFI N–S limit (hr) Aggregation by resampling Aggregation by moving average Cut-off point (hr) 
Blind Beck 0.25 6.3 22.1 47% 19.9 2.5 h (10 time periods) 2.25 h (9 time periods) 3.8 
Baru 0.083 1.1 18.7 62% 3.4 0.9–1 h (11–12 time periods) 1 h (12 time periods) 1.7 
            Time resolution estimates
 
  
Catchment Sampling frequency (hr) TCq (hr) TCs (hr) SFI N–S limit (hr) Aggregation by resampling Aggregation by moving average Cut-off point (hr) 
Blind Beck 0.25 6.3 22.1 47% 19.9 2.5 h (10 time periods) 2.25 h (9 time periods) 3.8 
Baru 0.083 1.1 18.7 62% 3.4 0.9–1 h (11–12 time periods) 1 h (12 time periods) 1.7 
Figure 5

Comparison of aggregated sequence to the IER sequence for (a) Blind Beck (sampling interval 15 min) and (b) Baru (sampling interval 5 min) at aggregations of 4, 8, 12 and 24 time periods (samples) illustrating how aggregation lowers the peak and spreads the volume of rainfall over a longer time period. The IER sequence is plotted for comparison.

Figure 5

Comparison of aggregated sequence to the IER sequence for (a) Blind Beck (sampling interval 15 min) and (b) Baru (sampling interval 5 min) at aggregations of 4, 8, 12 and 24 time periods (samples) illustrating how aggregation lowers the peak and spreads the volume of rainfall over a longer time period. The IER sequence is plotted for comparison.

Figure 6

Rt2 and R tend to a maximum value as aggregation increases for (a) Blind Beck and (b) Baru. The resolution of the IER is taken to be the point at which the maximum is reached or very little change is apparent. For Blind Beck, this value is reached at 10 periods for both Rt2 and R. The result for Baru is not quite as clear but can be estimated to be 10 periods from R and 11 or 12 from Rt2 although Rt2 continues to increase up to 24 time periods, perhaps due to higher variability of the rainfall.

Figure 6

Rt2 and R tend to a maximum value as aggregation increases for (a) Blind Beck and (b) Baru. The resolution of the IER is taken to be the point at which the maximum is reached or very little change is apparent. For Blind Beck, this value is reached at 10 periods for both Rt2 and R. The result for Baru is not quite as clear but can be estimated to be 10 periods from R and 11 or 12 from Rt2 although Rt2 continues to increase up to 24 time periods, perhaps due to higher variability of the rainfall.

Figure 7

A similar plot to Figure 6 with aggregation by moving average for (a) Blind Beck and (b) Baru. Rather than reaching an asymptotic level, the Rt2 and R values maximise at 9 time periods for Blind Beck and 12 time periods for Baru (determined graphically in Matlab). These values have been used as the estimates of the resolution of the IER and agree well with the estimates made by resampling. Convolution term in the caption is with reference to the method of calculating the moving average.

Figure 7

A similar plot to Figure 6 with aggregation by moving average for (a) Blind Beck and (b) Baru. Rather than reaching an asymptotic level, the Rt2 and R values maximise at 9 time periods for Blind Beck and 12 time periods for Baru (determined graphically in Matlab). These values have been used as the estimates of the resolution of the IER and agree well with the estimates made by resampling. Convolution term in the caption is with reference to the method of calculating the moving average.

Table 1 shows that the estimated resolution of the IER sequence for Blind Beck is around 9–10 time periods (i.e., 2.25–2.5 h) and for Baru it is 11–12 time periods (i.e., 55 min–1 hr). Both estimates are within the N–S safe sampling limit and below the fast time constant for both catchments indicating that even though resolution has been lost – the regularisation trade-off for numerical stability – the dominant mode of the rainfall–streamflow dynamics has been captured. Table 2 shows that the estimated resolution of the IER for both catchments is well within the Nyquist limit and, while the Blind Beck resolution is within the safe limits suggested by Ljung (1999) and Young (2010), the estimated resolution for Baru is close to the fast TC and outside the suggested limits. The estimates of resolution of the inferred sequence made from the aggregation plots are not always well defined and may be dependent on the length of record which will affect the number of aggregation periods that may be meaningfully calculated given the finite length of the data series. A better means of estimation of resolution may be achieved by examining the frequency spectra of the rainfall and streamflow sequences.

Table 2

The estimated resolution of the IER compared to the TCq, Nyquist limit and safe sampling limits suggested by Ljung (1999) and Young (2010) 

Catchment TCq (hours) Nyquist limit (hr) Ljung interval (hr) Young interval (hr) Estimated resolution (hr) 
Blind Beck 6.3 19.9 3.98 3.32 2.25–2.5 
Baru 1.1 3.4 0.68 0.57 0.91–1.0 
Catchment TCq (hours) Nyquist limit (hr) Ljung interval (hr) Young interval (hr) Estimated resolution (hr) 
Blind Beck 6.3 19.9 3.98 3.32 2.25–2.5 
Baru 1.1 3.4 0.68 0.57 0.91–1.0 

In Figure 8, the amplitude spectra of inferred effective and OER are very close (overlapping when smoothed) within a broad range of frequencies. The cut-off frequency, where the difference between the smoothed ER and IER spectra is approximately −6 dB, provides a frequency domain estimate of the resolution. The cut-off period for Blind Beck is 3.8 h and for Baru is 1.7 h. For frequencies above this value, a very strong low-pass filtering effect shown is by the rapid decrease in the IR spectrum. The frequency range beyond the cut-off point, shaded in Figure 8, carries a very small proportion of the power of the signal and can be considered non-significant.

Figure 8

Periodograms for (a) Blind Beck and (b) Baru showing the frequency structure of the effective rain, inferred effective rain and streamflow sequences. The grey area shows frequencies beyond the 6 dB difference between smoothed ER and IER spectral point. Both catchments show a similarity in the frequency spectra of effective and IER within the catchment system. The IER spectrum is very close to that of the actual effective rainfall within a wide range of frequencies mostly covering those corresponding to the catchments' time constants. There is also a strong low-pass filtering effect cutting off high frequencies with low amplitudes instead of boosting this high frequency noise.

Figure 8

Periodograms for (a) Blind Beck and (b) Baru showing the frequency structure of the effective rain, inferred effective rain and streamflow sequences. The grey area shows frequencies beyond the 6 dB difference between smoothed ER and IER spectral point. Both catchments show a similarity in the frequency spectra of effective and IER within the catchment system. The IER spectrum is very close to that of the actual effective rainfall within a wide range of frequencies mostly covering those corresponding to the catchments' time constants. There is also a strong low-pass filtering effect cutting off high frequencies with low amplitudes instead of boosting this high frequency noise.

Table 1 lists the time constants, SFIs and cut-off points for both catchments. The cut-off point for Blind Beck (3.8 h) is outside the range of the catchment time constants (6.3–22.1 h), probably reflecting the frontal rainfall regime, which is relatively uniform in time and space. Flow along both pathways is almost evenly split indicating that they are both important in terms of flow generation. On the other hand, Baru's TCq (1.1 h) is beyond the cut-off point (1.7 h) in the area where the spectra contain little power or information indicating why the catchment's variable, high intensity, high frequency, highly localised convective rainfall may not be easy to estimate. It is worth noting that the forward rainfall–discharge model does not fit the Baru catchment (88%) characterised by its highly variable, both spatially and temporally rainfall, as well as the Blind Beck catchment (98%) with its relatively uniform predominantly frontal rainfall (Kretzschmar et al. 2014).

The processes and characteristics limiting the IER accuracy include the slow components of the catchment dynamics and the rainfall regime. These can be seen as the ‘usual suspects’ affecting the inversion process. The general goodness-of-fit of the initial catchment model (rainfall–streamflow) appears to be a factor as well (see Figure 9), indicating that the IER estimation method presented here can be used to assess the quality of available data and the degree to which the data characterise the catchment. Further work is required using a range of catchment and rainfall regimes to confirm these results and explain them in terms of rainfall and catchment characteristics, as well as investigating spatial relationships. The latter will be evaluated using catchment data with multiple rain-gauges, and are the subject of the forthcoming work. Rainfall is the key driver of streamflow with the pattern varying from event to event, however, the underlying catchment characteristics, for example, soil, geology, topography, may modify this. A combination of inversion and spectral analysis may provide a method for untangling the effects of catchment characteristics and rainfall regime on streamflow generation and has the potential for characterising the effects of future changes in catchment and/or rainfall characteristics due to, for example, climate change.

Figure 9

Comparison of observed discharge and discharge generated from the IER for (a) Blind Beck and (b) Baru. The flow sequences match almost perfectly in the case of Blind Beck and very closely in the case of Baru where the peak flows are underestimated. Note that the forward (rainfall–discharge) model fit for Blind Beck (98%) is better than for Baru (88%) (Kretzschmar et al. 2014).

Figure 9

Comparison of observed discharge and discharge generated from the IER for (a) Blind Beck and (b) Baru. The flow sequences match almost perfectly in the case of Blind Beck and very closely in the case of Baru where the peak flows are underestimated. Note that the forward (rainfall–discharge) model fit for Blind Beck (98%) is better than for Baru (88%) (Kretzschmar et al. 2014).

CONCLUSIONS

A combination of time and frequency domain techniques has been used to show that the IER time-series generated by the RegDer inversion method does, indeed, approximate the direct inverse of a transfer function to a high degree of accuracy within the frequency range which includes the dominant modes of the rainfall–streamflow dynamics. The direct inverse exaggerates low-amplitude high frequency noise, which is filtered out by the regularisation process involved in the RegDer method. The smoothing of the signal resulting from regularisation is quantified in the time domain by comparison with aggregated observed input data using standard model fit measures – coefficient of determination, Rt2, and correlation coefficient, R – and analysed as a low-pass filtering process in the frequency domain. The smoothing effect is minimised within the constraints of the available data and catchment dynamics, through optimisation of the regularisation constants (NVRs) to obtain the best fit of the inversion process where both rainfall and discharge data are available.

ACKNOWLEDGEMENTS

The authors would like to thank Mary Ockenden for the collection and quality assurance of the period of rainfall and streamflow for the Blind Beck catchment (NERC grant number NER/S/A/2006/14326), and also Jamal Mohd Hanapi and Johnny Larenus for the collection of the period of rainfall and streamflow utilised for the Baru catchment and to Paul McKenna for its quality assurance (NERC grant number GR3/9439). This work has been partly supported by the Natural Environment Research Council (Consortium on Risk in the Environment: Diagnostics, Integration, Benchmarking, Learning and Elicitation (CREDIBLE)) grant number: NE/J017299/1.

REFERENCES

REFERENCES
Andrews
F.
Croke
B.
Jeanes
K.
2010
Robust Estimation of the Total Unit Hydrograph
.
In: 2010 International Congress on Environmental Modelling and Software Modelling for Environment's Sake. Ottawa
,
Canada
.
Beven
K. J.
2012
Rainfall–Runoff Modelling: The Primer, 2nd edn.
John Wiley and Sons
,
Chichester, UK
.
Bloomfield
P.
1976
Fourier Analysis of Time Series: An Introduction
.
John Wiley & Sons
,
New York, USA
.
Blöschl
G.
Sivapalan
M.
Wagener
T.
Viglione
A.
Savenije
H.
(eds)
2013
Runoff Prediction in Ungauged Basins: Synthesis across Processes, Places and Scales
.
Cambridge University Press
,
New York, USA
.
Boorman
D.
1989
A new approach to Unit Hydrograph Modelling
.
PhD thesis
,
Lancaster University
,
Lancaster, UK
.
Brocca
L.
Moramarco
T.
Melone
F.
Wagner
W.
2013
A new method for rainfall estimation through soil moisture observations
.
Geophys. Res. Lett.
40
(
5
),
853
858
.
Brocca
L.
Ciabatta
L.
Massari
C.
Moramarco
T.
Hahn
S.
Hasenauer
S.
Kidd
R.
Dorigo
W.
Wagner
W.
Levizzani
V.
2014
Soil as a natural rain gauge: estimating global rainfall from satellite soil moisture data
.
J. Geophys. Res. Atmos.
119
(
9
),
5128
5141
.
Chappell
N.
Franks
S.
Larenus
J.
1998
Multi-scale permeability estimation for a tropical catchment
.
Hydrol. Process.
12
(
9
),
1507
1523
.
Chappell
N. A.
McKenna
P.
Bidin
K.
Douglas
I.
Walsh
R. P. D.
1999
Parsimonious modelling of water and suspended-sediment flux from nested catchments affected by selective tropical forestry
.
Philos. Trans. Royal Soc. Lond. Ser. B.
354
,
1831
1846
.
Chappell
N. A.
Tych
W.
Chotai
A.
Bidin
K.
Sinun
W.
Chiew
T. H.
2006
Barumodel: combined data based mechanistic models of runoff response in a managed rainforest catchment
.
For. Ecol. Manage.
224
(
1
),
58
80
.
Chappell
N. A.
Sherlock
M.
Bidin
K.
Macdonald
R.
Najman
Y.
Davies
G.
Sawada
H.
Araki
M.
LaFrankie
J. V.
Shimizu
A.
2007
Runoff processes in Southeast Asia: Role of soil, regolith and rock type
. In:
Forest Environments in the Mekong River Basin
(
Sawada
H.
Araki
M.
Chappell
N. A.
LaFrankie
J. V.
Shimizu
A.
, eds),
Springer
,
Tokyo
,
Japan
, pp.
3
23
.
Croke
B.
2010
Exploring changes in catchment response characteristics: application of a generic filter for estimating the effective rainfall and unit hydrograph from an observed streamflow time series
. In:
Proceedings, British Hydrological Society International Conference: Role of Hydrology in Managing Consequences of a Changing Global Environment
,
July 19–23, 2010
,
Newcastle University
,
Newcastle upon Tyne
,
UK
, pp.
494
499
.
Cuchi
J.
Chinarro
D.
Villarroel
J.
Antonio Cuchi
D.
Luis Villarroel
J.
2014
Linear system techniques applied to the Fuenmayor Karst Spring, Huesca (Spain)
.
Environ. Earth Sci.
71
(
3
),
1049
1060
.
Holman-Dodds
J. K.
Bradley
A. A.
Sturdevant-Rees
P. L.
1999
Effect of temporal sampling of precipitation on hydrologic model calibration
.
J. Geophys. Res.
104
(
D16
),
19645
19654
.
Krajewski
W. F.
Lakshmi
V.
Georgakakos
K. P.
Jain
S. C.
1991
A Monte Carlo study of rainfall sampling effect on a distributed catchment model
.
Water Resour. Res.
27
(
1
),
119
128
.
Kretzschmar
A.
Tych
W.
Chappell
N. A.
2014
Reversing hydrology: estimation of sub-hourly rainfall time-series from streamflow
.
Environ. Modell. Softw.
60
,
290
301
.
Krier
R.
Matgen
P.
Goergen
K.
Pfister
L.
Hoffmann
L.
Kirchner
J. W.
Uhlenbrook
S.
Savenije
H. H. G.
2012
Inferring catchment precipitation by doing hydrology backward: a test in 24 small and mesoscale catchments in Luxembourg
.
Water Resour. Res.
48
,
W10525
.
Littlewood
I. G.
2007
Rainfall–streamflow models for ungauged basins: uncertainty due to modelling time step
. In:
Uncertainties in the ‘Monitoring–Conceptualisation– Modelling’ Sequence of Catchment Research (Proc. Eleventh Biennial Conference of the Euromediterranean Network of Experimental and Representative Basins)
(L. Pfister & L. Hoffmann, eds)
.
UNESCO Tech. Doc. in Hydrology Series 81
,
Paris
,
France
, pp.
149
155
.
Littlewood
I. G.
Young
P. C.
Croke
B. F. W.
2010
Preliminary comparison of two methods for identifying rainfall–streamflow model parameters insensitive to data time-step: the Wye at Cefn Brwyn, Plynlimon, Wales
. In:
Proceedings of the Third International Symposium (British Hydrological Society)
,
19–23 July 2010
,
Newcastle University
,
UK
.
Ljung
L.
1999
System Identification. Theory for the User
,
2nd edn
.
Prentice Hall
,
Upper Saddle River, NJ
,
USA
.
Mayes
W. M.
Walsh
C. L.
Bathurst
J. C.
Kilsby
C. G.
Quinn
R. F.
Wilkinson
M. E.
Daugherty
A. J.
Connell
P. E.
2006
Monitoring a flood event in a densely instrumented catchment, the Upper Eden, Cumbria, UK
.
Water Environ. J.
20
(
4
),
217
226
.
Reynard
N. S.
Stewart
E. J.
1993
The derivation of design rainfall profiles for upland areas of the UK
.
Meteorol. Mag.
122
,
116
123
.
Taylor
C. J.
Pedregal
D. J.
Young
P. C.
Tych
W.
2007
Environmental time series analysis and forecasting with the Captain toolbox
.
Environ. Model. Softw.
22
(
6
),
797
814
.
Teuling
A. J.
Lehner
I.
Kirchner
J. W.
Seneviratne
S. I.
2010
Catchments as simple dynamical systems: experience from a Swiss pre-alpine catchment
.
Water Resour. Res.
46
(
10
),
W10502
.
Walsh
R. P. D.
Bidin
K.
Blake
W. H.
Chappell
N. A.
Clarke
M. A.
Douglas
I.
Ghazali
R.
Sayer
A. M.
Suhaimi
J.
Tych
W.
Annammala
K. V.
2011
Long-term responses of rainforest erosional systems at different spatial scales to selective logging and climatic change
.
Phil. Trans. R. Soc. B
366
,
3340
3353
.
Wickert
M.
2013
Signals and Systems for Dummies
,
1st edn
.
John Wiley & Sons
,
Hoboken, NJ
,
USA
.
Young
P. C.
2010
The estimation of continuous-time rainfall–flow models for flood risk management
. In:
Proceedings of the Third International Symposium (British Hydrological Society)
,
19–23 July 2010
,
Newcastle University
,
UK
.
Young
P. C.
Romanowicz
R. J.
2004
PUB and data-based mechanistic modelling: the importance of parsimonious continuous-time models. In:
Proceedings of the iEMSs 2004 International Congress: complexity and integrated resources management
.
International Environmental Modelling and Software Society
,
Osnabruech, Germany
, pp.
214
224
.
Young
P. C.
Sumisławska
M. A.
2012
Control systems approach to input estimation with hydrological applications
. In:
16th IFAC Symposium on System Identification
.
11–13 July
,
Brussels
,
Belgium
, pp.
1043
1048
.
Young
P. C.
Pedregal
D. J.
Tych
W.
1999
Dynamic harmonic regression
.
J. Forecast.
18
(
6
),
369
394
.