# Imputation of Missing Streamflow Data at Multiple Gauging Stations in Benin
Republic

Rendani Mbuvha,12 Peniel Julien Yise Adounkpe, 3 Wilson Tsakane Mongwe,4
Mandela Coovi Mahuwetin Houngnibo, 5 Nathaniel Newlands,6 Tshilidzi Marwala 4

###### Abstract

Streamflow observation data is vital for flood monitoring, agricultural, and
settlement planning. However, such streamflow data are commonly plagued with
missing observations due to various causes such as harsh environmental
conditions and constrained operational resources. This problem is often more
pervasive in under-resourced areas such as Sub-Saharan Africa. In this work,
we reconstruct streamflow time series data through bias correction of the
GEOGloWS ECMWF streamflow service (GESS) forecasts at ten river gauging
stations in Benin Republic. We perform bias correction by fitting Quantile
Mapping, Gaussian Process, and Elastic Net regression in a constrained
training period. We show by simulating missingness in a testing period that
GESS forecasts have a significant bias that results in low predictive skill
over the ten Beninese stations. Our findings suggest that overall bias
correction by Elastic Net and Gaussian Process regression achieves superior
skill relative to traditional imputation by Random Forest, k-Nearest
Neighbour, and GESS lookup. The findings of this work provide a basis for
integrating global GESS streamflow data into operational early-warning
decision-making systems (e.g., flood alert) in countries vulnerable to drought
and flooding due to extreme weather events.

## Introduction

Discharge predictions are critical for decision-making in many areas of
climate change adaptation, including flood and drought prevention,
agricultural planning, and hydroelectric power system operations (Kratzert et
al. 2019). Predictions can be generated by physical and statistical models of
the hydrological process within the river basin (Adounkpe et al. 2021).
Regardless of the approach, observation data is required to sufficiently
calibrate the localised dynamics within catchments. Model-based predictions
are even more essential in vulnerable and under-developed regions that are
less resilient to extreme weather events such as drought and flooding.
However, localised observations for hydrological models are not always
available due to physical sensor failure events (Hamzah et al. 2020). Physical
sensor systems often endure harsh environmental conditions, destruction and
power outages (Tencaliec et al. 2015). These phenomena result in prolonged
periods of missing or invalid observation readings. In developing regions such
as Benin Republic, the effect of such missing data is further complemented by
the already sparse distribution of stream gauging stations (Harrigan et al.
2020a). Fitting statistical models based on incomplete observation data
results in inaccurate models that do not capture the full dynamics of the
underlying hydrological processes (Hamzah et al. 2020). A more complete set of
streamflow data allows for greater confidence in the outputs of hydrological
models (Pigott 2001). Impuation of missing data therefore leads to much needed
improved flood and drought predictions in a climate change context where these
events are expected to increase in intensity and frequency (Kalantari et al.
2018). Thus missing data imputation has become a necessary first step in data
preprocessing tasks (Pigott 2001). Moreover, missing observations and the
estimation of streamflow at ungauged locations are two of the most pervasive
problems in the field of hydrology (Hamzah et al. 2020). Since the area of
missing data imputation commands a relatively vast amount of statistical
literature, a wide variety of techniques have been developed and applied with
different data requirements and accuracy (Little and Rubin 2019).

A simple approach to handling missing data is discarding records where any
missingness arises. Gill et al. (2007) show that while deletion of missing
data is common practice in hydrology, it is not necessarily optimal as
significant amounts of predictive information can be lost. On the other hand,
single imputation methods attempt to replace missing values one feature at a
time with a pre-specified summary statistic (e.g. mean or median) of the
complete data (Baraldi and Enders 2010). Regression-based methods such as
imputation by Random Forest, k-Nearest Neighbours and Multilayer perceptrons
improve upon single imputation by exploiting correlations between variables in
complete data (Gao 2017; Pantanowitz and Marwala 2009). Both regression and
single imputation can lead to biases in the resultant analysis when the
assumption that data are missing at random is violated (Gao 2017). Multiple
imputation methods such as Multivariate Imputation by Chained Equations (MICE)
reduce this bias by replacing missing data points with predictions from an
ensemble of regressors obtained from multiple subsets of the complete data
(Van Buuren 2018).

Each of the imputation techniques above is used extensively in hydrology.
Hamzah et al. (2021) compares numerous regression methods to reconstruct
streamflow data in the Malaysian Langat River Basin. Norazizi and Deni (2019)
employs MICE and expectation maximisation (EM) for the imputation of rainfall
data. Adeloye and Rustum (2012) and Mwale, Adeloye, and Rustum (2012) propose
the use of self-organising maps (SOMs) for infilling streamflow data at
inadequately gauged basins in South West Nigeria and the Shine river basin in
Malawi. In this work, we propose further augmenting regression imputation with
data from Global streamflow forecasting systems.

Global streamflow prediction systems present a potential solution to missing
data and the sparse distribution of stream gauge stations in developing
countries. The Global Flood Awareness System (GLoFAS) is one such system
jointly developed by the European Commission’s Joint Research Center (JRC) and
the European Centre for Medium-Range Weather Forecasting (ECMWF) (Harrigan et
al. 2020b). GloFAS transforms numerical weather predictions to surface runoff
through the hydrology tiled ECMWF scheme for surface exchanges over land
(HTESSEL) model (Sanchez Lozano et al. 2021). These surface runoff estimates
are then routed to stream networks using the LISFOOD river routing scheme
(Snow 2015; Harrigan et al. 2020b). A limitation of the GloFAS system is that
streamflow estimates are only available at a low resolution of a 0.1-degree
grid cell (10 km2); this limits the applicability of GloFAS to mainly large
watersheds (Snow 2015; Souffront Alcantara et al. 2019) .

The Group on Earth Observations Global Water and Sustainability Initiative
(GEOGloWS) attempts to address the coarse resolution of GloFAS with the
GEOGloWS ECMWF streamflow service (GESS) (Ashby et al. 2021). The GESS extends
GloFAS to local sites using an updated digital elevation model (DEM) that
allows for increased resolutions in areas of the world with narrow watersheds
(Ashby et al. 2021). Similar to GloFAS, runoff predictions from the HTESSEL
land surface model are further downscaled and routed through a high-resolution
stream network using the Muskingum algorithm of the routing application for
parallel-discharge computation (RAPID) (Gill 1978).

While GloFAS and GESS provide much-needed streamflow data, localised
calibration is still a significant challenge. A global calibration of these
systems is biased primarily to developed regions with dense and complete
observations data (Ashby et al. 2021). Thus, the utility of streamflow
forecasts from these systems can be significantly enhanced by bias and
variance correction to minimise systemic biases between model output and
gauging station observations. Statistical Learning methods are broadly popular
when learning transfer functions between global model forecasts and local
observations (Piani et al. 2010; Rischard, Pillai, and McKinnon 2018; Wang,
Zhang, and Villarini 2021). Various statistical learning methods have been put
forward for bias correction tasks, including penalised regression, Quantile
Mapping, Artificial Neural Networks and Gaussian Processes (Xu and Liang 2021;
Pastén-Zapata et al. 2020; Hunt et al. 2022).

This paper presents and evaluates the bias correction of GESS river discharge
estimates using statistical learning methods on ten hydrological gauging
stations in Benin Republic. The proposed bias correction methods can be
utilised to impute missing observation data with higher accuracy than
standalone GESS and are superior to complete data-based Random Forest and
k-Nearest Neighbour imputation.

Figure 1: Benin’s catchments, rivers and hydrological stations missing data
rate.

## Methods

Figure 2: Missing data periods at each hydrological gauging station.

### Study Area

For the purposes of this study, ten hydrological stations located on the
outlets of Benin’s major river basins were selected: Athiémè, Beterou, Bonou,
Domè, Kaboua, Kompongou, Koubéri, Lanta, Porga and Yakin. Benin’s basins are
not limited to only Benin’s administrative boundary but extend to its
neighbouring countries: Togo, Burkina Faso and Nigeria. The Kouffo, Mono,
Okpara, Ouémé and Zou rivers originate from central Benin and Togo and
generally flow into the Atlantic Ocean. The Alibori, Mékrou and Sota rivers
also originate from the centre of Benin but have their outlet up north in the
Niger river. The Pendjari river (sometimes referred to as the Oti river)
originates from Northern Benin and flows into the Volta river in Togo and
Ghana.

### Data

In-situ hydrological data (daily river discharge record at ten locations from
1980 to 2021) was acquired from Benin’s hydrological service, Direction
Générale de l’Eau (DG-Eau). The percentage of missing data for the selected
hydrological stations is 27.9%, with the station of Athiémè having the least
and the Kompongou station having the most. Figure 2 shows in detail the
periods of missing data for each station. On the other hand, the GESS
forecasts present no missing data and provide daily river discharge data
ranging from 1979 to the present at targeted river sections (Ashby et al.
2021).

Figure 3: Plot of Mean KGE, NSE and RMSE measures from each imputation method
at varying levels of missingness across the ten stations.

### Imputation Methods

We investigate streamflow imputation by bias-correcting GESS forecasts. We
propose an imputation scheme where we replace a missing in-situ observation at
time $t$ with:

$\tilde{\boldsymbol{x}_{t}}=h(\boldsymbol{x}_{t}^{\prime})$ (1)

where $\boldsymbol{x}^{\prime}$ the GESS forecast of streamflow at time $t$.
$h$ is a regressor trained on periods where $\boldsymbol{x}$ and
$\boldsymbol{x}^{\prime}$ are both available. We investigate settings where
$h$ takes various forms, including; a simple lookup from GESS, a quantile
mapping, Elastic Net and Gaussian Process.

#### Gaussian processes (GPs)

are nonparametric probabilistic models, that are well known for regression
problems (Cohen et al. 2020). A GP is a collection of random variables, where
each of its finite subsets is jointly Gaussian (Rasmussen and Williams 2006).
GPs are sufficiently defined by a mean function $\mu(\boldsymbol{x}^{\prime})$
and a kernel function
$k(\boldsymbol{x}_{i}^{\prime},\boldsymbol{x}_{j}^{\prime})$.

Given a training dataset
$\\{\mathbf{x}^{\prime}_{t},\mathbf{x}_{t}\\}^{t}_{t=1}$ with $t$ noisy
observations $\mathbf{x}_{t}=h(\boldsymbol{x}^{\prime}_{t})+\epsilon$, where
$\epsilon\sim\mathcal{N}(0,\sigma_{X}^{2})$. The GP prior on $h$ is such that
$h(\boldsymbol{x}^{\prime})\sim\mathcal{N}(\boldsymbol{\mu}_{x^{\prime}},\boldsymbol{K}_{x^{\prime}}+\sigma_{X}^{2}\boldsymbol{I})$.
The mean and variance of the Gaussian posterior predictive distribution of the
function $h(\boldsymbol{x^{\prime}}_{*})$ at a test point
$\boldsymbol{x}_{*}^{\prime}$, are given by

$\displaystyle\mathbb{E}(h_{*})$
$\displaystyle=\mathbf{\mu_{x^{\prime}_{*}}}+\boldsymbol{k}_{*}^{T}(\boldsymbol{K}_{x^{\prime}}+\sigma_{X^{\prime}}^{2}\mathbf{I})^{-1}(\mathbf{X}-\mathbf{\mu_{x^{\prime}})},$
(2) $\displaystyle\text{var}(h_{*})$
$\displaystyle=\boldsymbol{k}_{**}-\mathbf{k}_{*}^{T}(\boldsymbol{K}_{x^{\prime}}+\sigma_{X^{\prime}}^{2}\mathbf{I})^{-1}\mathbf{k}_{*}$
(3)

where $\boldsymbol{k}_{**}=k(\mathbf{x^{\prime}}_{*},\mathbf{x^{\prime}}_{*})$
and $\boldsymbol{k}_{*}=k(\mathbf{x^{\prime}},\mathbf{x^{\prime}}_{*})$. The
kernel hyperparameters $\boldsymbol{\theta}$ and the noise parameter
$\sigma_{\boldsymbol{X}}$ are obtained by maximising the log-marginal
likelihood

$\log
p(\boldsymbol{X}|\boldsymbol{X}^{\prime})=\log\mathcal{N}(\boldsymbol{\mu}_{x^{\prime}},\boldsymbol{K}_{x^{\prime}x^{\prime}}+\sigma_{\boldsymbol{X}}^{2}\mathbf{I}\big{)}.$
(4)

In this work, we use a multi-output formulation of the GP where all ten
stations share a squared exponential kernel, thus allowing for implicit
inference of connectivity relationships between stations.

#### Elastic Net

The elastic net is a penalised linear regression model that employs weighted
$l_{1}$ and $l_{2}$ norm regularisation terms. The regularisation terms serve
as conservative priors (bias towards zero) on the coefficients that prevent
over-fitting to training data. Similar to the GP, we utilised a multi-output
formulation of the Elastic Net that allows information sharing between
stations. In our imputation formulation;

$\tilde{\boldsymbol{x}}_{t}=h(\boldsymbol{x}^{\prime}_{t})=\boldsymbol{X}^{\prime}\boldsymbol{\beta}+\epsilon$
(5)

where $\boldsymbol{\beta}$ is a matrix of coefficients on the GESS predictions
$\boldsymbol{x}^{\prime}_{t}$. The matrix of coefficients is obtained by
optimising the posterior :

$\|\boldsymbol{X}-\boldsymbol{X}^{\prime}\boldsymbol{\beta}\|^{2}+\lambda_{2}\|\boldsymbol{\beta}\|^{2}+\lambda_{1}\|\boldsymbol{\beta}\|_{1}$
(6)

with the hyperparametes $\lambda_{1}$ and $\lambda_{2}$ obtained by cross
validation with the constraint $\lambda_{1}+\lambda_{2}=1$. Given limited in-
situ data, the elastic net has the advantage of limiting over-fitting the bias
correction relative to a regularised regression.

#### Quantile Mapping

Quantile Mapping (QM) is a well-known method for bias correction of physical
model output (Maraun 2013; Ringard, Seyler, and Linguet 2017). In QM, we seek
to learn a mapping between the empirical cumulative distribution functions
(CDF) of the in-situ and GESS prediction. In this case, the bias correction
transfer function is:

$\tilde{\boldsymbol{x}}_{t}=h(\boldsymbol{x}^{\prime}_{t})=F_{X}^{-1}F_{X^{\prime}}(x^{\prime}_{t})$
(7)

Where $F_{X}^{-1}$ is the inverse empirical CDF of the in-situ data and
$F_{X^{\prime}}$ is the empirical CDF of the GESS data. Both empirical CDFs
are obtained during a specified training period.

####  Baselines

We compare imputation by bias correction using the abovementioned methods to
traditional regression imputation methods that rely only on complete in-situ
data. We consider regression-based imputation by Random Forest (RF)
(Pantanowitz and Marwala 2009) and k-nearest neighbours Hamzah et al. (2021)
as our baselines.

### Performance Evaluation

We measure the quality of the imputation using Kling-Gupta Efficiency (KGE)
(Gupta et al. 2009), Nash-Sutcliffe Efficiency (NSE) (Nash and Sutcliffe 1970)
and Root Mean Square Error (RMSE) (Chai and Draxler 2014).

The KGE metric provides a diagnostically informative decomposition of the NSE
and RMSE. It facilitates the analysis of the relative importance of its
different components (correlation, bias and variability) in the context of
hydrological modelling. The NSE measures the ability to predict variables
different from the mean and gives the proportion of the initial variance
accounted for by the model. The RMSE is frequently used to evaluate how
closely the predicted values match the observed values based on the relative
range of the data.

$\displaystyle
KGE=1-\sqrt{(r-1)+\left(\dfrac{\bar{\hat{x_{t}}}}{\bar{x_{t}}}-1\right)^{2}+\left(\dfrac{\sigma_{\hat{x_{t}}}}{\sigma_{x_{t}}}-1\right)^{2}}$
(8)

$\displaystyle
NSE=1-\dfrac{\sum\limits_{i=1}^{n}(x_{t}-\hat{x_{t}})^{2}}{\sum\limits_{i=1}^{n}(x_{t}-\bar{x_{t}})^{2}}$
(9)

$\displaystyle
RMSE=\sqrt{\dfrac{1}{n}\sum_{i=1}^{n}\left(x_{t}-\hat{x_{t}}\right)^{2}}$ (10)

In the above equations, $x_{t}$ represents the observed discharge at a time
$t$, $\hat{x_{t}}$ represents the imputed discharge at a time $t$,
$\bar{x_{t}}$ is the mean of the observed discharge, $\bar{\hat{x_{t}}}$ is
the mean of the imputed discharge, $\sigma_{\hat{x_{t}}}$ is the variance of
the simulated discharge, $\sigma_{\hat{x_{t}}}$ is the variance of the
observed discharge, $r$ is the Pearson correlation coefficient, and $n$ is the
number of observations. The KGE values range from $-\infty$ to 1, with values
close to 1 indicating good agreement. Similarly to the KGE, the NSE values
also range from $-\infty$ to 1, with the perfect model having the value 1. The
RMSE values range from 0 to $\infty$, with 0 implying the model has a perfect
fit.

## Results and Discussion

We partition the data such that 60% of the period (1980/01-2004/11) is for
training and 40% is for (2004/11-2021/06) testing. We evaluate the performance
of the imputation methods by randomly simulating missingness on the complete
testing data at rates of 5%, 10%, 20%, 30% and 50%.

Figure 3 shows the mean KGE, NSE and RMSE of the imputation methods across
gauging stations and levels of missingness. It can be seen that a simple
lookup of the GESS predictions yields the worst imputation performance. The
low KGE and high RMSE values detailed in Table 1 suggest that the GESS
predictions have significant bias.

The hypothesis around significant bias in the GESS predictions is reinforced
by the outperformance of the GP and Elastic Net bias correction imputation.
Importantly the imputation by bias correction yields better performance than
simply employing RF, KNN and MLP regressors on the complete in-situ data. This
can be reasonably expected as there are limited periods where the data is
jointly complete for training the multiple imputation regressors. Thus the
GESS data provides the necessary signals regarding discharge seasonality and
hydrologic connectivity to enhance imputation quality significantly.

Method | NSE | RMSE
---|---|---
% Missing Data
| 5% | 10% | 20% | 30% | 50% | 5% | 10% | 20% | 30% | 50%
KNN | 0.950 | 0.892 | 0.739 | 0.596 | 0.389 | 24.694 | 36.802 | 57.686 | 71.743 | 88.948
RF | 0.947 | 0.888 | 0.736 | 0.593 | 0.372 | 25.222 | 37.476 | 58.400 | 72.625 | 90.265
Gaussian Process | 0.973 | 0.938 | 0.840 | 0.751 | 0.623 | 17.998 | 28.308 | 45.950 | 57.364 | 71.025
GESS | 0.799 | 0.655 | 0.447 | 0.121 | -0.378 | 38.219 | 59.650 | 79.612 | 105.145 | 129.798
Elastic Net | 0.973 | 0.948 | 0.872 | 0.804 | 0.706 | 17.783 | 24.484 | 38.719 | 48.309 | 59.559
Quantile Mapping | 0.966 | 0.938 | 0.859 | 0.787 | 0.671 | 19.557 | 26.900 | 40.323 | 49.993 | 61.832

Table 1: Mean KGE, NSE and RMSE of each imputation method at varying levels of
missingness across the ten stations. Figure 4: The infilled streamflow time
series at the Koubéri station using the GESS lookup Bias Correction.

Table 2 shows the KGE values at each of the stations at a level 20%
missingness. It can be seen that GESS yields the best performance at Yankin
and Kompongou stations. Yankin and Kompongou had the highest prevalence of
missing data at 47.67% and 59.48%, respectively. A possible reason for this
result is that, given the elevated levels of missingness at these stations,
the complete data could not sufficiently calibrate the bias correction. GP and
Elastic Net imputation methods provide the best performance in eight of the
ten stations with the highest KGE values obtained at stations with highly
complete data in the Athiémè and Bonou stations. QM as bias correction yields
KGEs superior to multiple imputations by k-NN and RF, suggesting that bias
correction would remain the most efficient approach when constrained by
computational resources.

| Stations
---|---
Method | Yankin | Lanta | Kompongou | Athiémè | Kaboua | Beterou | Bonou | Porga | Koubéri | Domè
KNN | 0.786 | 0.728 | 0.725 | 0.864 | 0.772 | 0.847 | 0.846 | 0.851 | 0.784 | 0.797
RF | 0.821 | 0.784 | 0.780 | 0.859 | 0.841 | 0.807 | 0.832 | 0.839 | 0.626 | 0.822
Gaussian Process | 0.863 | 0.857 | 0.829 | 0.919 | 0.833 | 0.831 | 0.877 | 0.913 | 0.914 | 0.894
GESS | 0.896 | 0.651 | 0.905 | 0.643 | 0.870 | 0.865 | 0.608 | 0.853 | 0.443 | 0.448
Elastic Net | 0.848 | 0.872 | 0.809 | 0.930 | 0.881 | 0.881 | 0.908 | 0.894 | 0.901 | 0.905
Quantile Mapping | 0.831 | 0.871 | 0.823 | 0.936 | 0.841 | 0.838 | 0.906 | 0.872 | 0.865 | 0.912

Table 2: KGE of each imputation method at 20% missingness at each respective
gauging station

Figure 4 depicts the embedded bias in the GESS forecasts at the Koubéri
station. Figure 4 shows the resultant time series by looking up the missing
data from the GESS. The positive bias in the GESS manifest in that GESS
lookups lead to new discharge extremes that are five times the extremes
recorded in the in-situ observations. In operational settings, such
significant positive biases can lead to false flood alerts and sub-optimal use
of already limited resources. The remedy to this extreme bias is seen in
Figure 5, where bias correction removes the significant positive bias while
retaining the seasonal periodicity in discharge.

Figure 5: The infilled streamflow time series at the Koubéri station using the
Elastic Net Bias Correction.

## Conclusion

We evaluated state-of-the-art imputation methods for streamflow prediction,
demonstrating the need and importance for bias-correcting GEOGloWS ECMWF
streamflow service (GESS) streamflow forecast information. Our findings show
that the bias introduced by GESS forecasts can be significant and result in
possible false flooding alerts. We minimize such bias using QM, GP, and
Elastic Net regressions trained in periods where in-situ observations are
available. The resultant imputation under simulated missingness shows that
GESS bias correction outperforms train k-NN and RF imputation trained on
available data alone.

Our findings enhance decision-making based on streamflow model-based
forecasts, reducing bias introduced by missing in-situ data. The reliance on
such forecasts is high, particularly in areas where data collection is too
costly, under-resourced, or not possible across Sub-Saharan Africa. The
proposed bias correction methodology integrates into existing operational GESS
forecasting and now enables real-time flooding and drought monitoring with
lead times of up to 15 days. Continued work aims to enhance and expand on the
operational implementation of our bias-correction-based methodology.
Specifically, we aim to interpolate and extrapolate the best-performing bias-
correction model to ungauged areas. This is anticipated to vastly increase the
utility of GESS to large areas across Sub-Saharan Africa and Central America
that are currently ungauged. The robustness of our proposed approach will also
be further tested at the seasonal scale to ensure that wet seasons and overall
seasonal dynamics are captured accurately.

## References

  * Adeloye and Rustum (2012) Adeloye, A. J.; and Rustum, R. 2012.  Self-organising map rainfall-runoff multivariate modelling for runoff reconstruction in inadequately gauged basins.  _Hydrology Research_ , 43(5): 603–617.
  * Adounkpe et al. (2021) Adounkpe, P. J. Y.; Alamou, E.; Diallo, B.; and Ali, A. 2021.  Predicting Discharge in Catchment Outlet Using Deep Learning: Case Study of the Ansongo-Niamey Basin.  In _NeurIPS 2021 Workshop on Tackling Climate Change with Machine Learning_.
  * Ashby et al. (2021) Ashby, K. R.; Hales, R. C.; Nelson, J.; Ames, D. P.; and Williams, G. P. 2021.  Hydroviewer: A Web Application to Localize Global Hydrologic Forecasts.  _Open Water Journal_ , 7(1): 9.
  * Baraldi and Enders (2010) Baraldi, A. N.; and Enders, C. K. 2010.  An introduction to modern missing data analyses.  _Journal of school psychology_ , 48(1): 5–37.
  * Chai and Draxler (2014) Chai, T.; and Draxler, R. R. 2014.  Root mean square error (RMSE) or mean absolute error (MAE)? – Arguments against avoiding RMSE in the literature.  _Geoscientific Model Development_ , 7(3): 1247–1250.
  * Cohen et al. (2020) Cohen, S.; Mbuvha, R.; Marwala, T.; and Deisenroth, M. 2020.  Healing products of Gaussian process experts.  In _International Conference on Machine Learning_ , 2068–2077. PMLR.
  * Gao (2017) Gao, Y. 2017.  _Dealing with missing data in hydrology: Data analysis of discharge and groundwater time-series in Northeast Germany_.  Ph.D. thesis, Freie Universität Berlin.
  * Gill (1978) Gill, M. A. 1978.  Flood routing by the Muskingum method.  _Journal of hydrology_ , 36(3-4): 353–363.
  * Gill et al. (2007) Gill, M. K.; Asefa, T.; Kaheil, Y.; and McKee, M. 2007.  Effect of missing data on performance of learning algorithms for hydrologic predictions: Implications to an imputation technique.  _Water resources research_ , 43(7).
  * Gupta et al. (2009) Gupta, H. V.; Kling, H.; Yilmaz, K. K.; and Martinez, G. F. 2009.  Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling.  _Journal of Hydrology_ , 377(1): 80–91.
  * Hamzah et al. (2021) Hamzah, F. B.; Hamzah, F. M.; Razali, S. M.; and Samad, H. 2021.  A comparison of multiple imputation methods for recovering missing data in hydrological studies.  _Civil Engineering Journal_ , 7(9): 1608–1619.
  * Hamzah et al. (2020) Hamzah, F. B.; Mohd Hamzah, F.; Mohd Razali, S. F.; Jaafar, O.; and Abdul Jamil, N. 2020.  Imputation methods for recovering streamflow observation: A methodological review.  _Cogent Environmental Science_ , 6(1): 1745133.
  * Harrigan et al. (2020a) Harrigan, S.; Zsoter, E.; Alfieri, L.; Prudhomme, C.; Salamon, P.; Wetterhall, F.; Barnard, C.; Cloke, H.; and Pappenberger, F. 2020a.  GloFAS-ERA5 operational global river discharge reanalysis 1979–present.  _Earth System Science Data_ , 12(3): 2043–2060.
  * Harrigan et al. (2020b) Harrigan, S.; Zsoter, E.; Alfieri, L.; Prudhomme, C.; Salamon, P.; Wetterhall, F.; Barnard, C.; Cloke, H.; and Pappenberger, F. 2020b.  GloFAS-ERA5 operational global river discharge reanalysis 1979–present.  _Earth System Science Data_ , 12(3): 2043–2060.
  * Hunt et al. (2022) Hunt, K. M. R.; Matthews, G. R.; Pappenberger, F.; and Prudhomme, C. 2022.  Using a long short-term memory (LSTM) neural network to boost river streamflow forecasts over the western United States.  _Hydrology and Earth System Sciences Discussions_ , 2022: 1–30.
  * Kalantari et al. (2018) Kalantari, Z.; Ferreira, C. S. S.; Keesstra, S.; and Destouni, G. 2018.  Nature-based solutions for flood-drought risk mitigation in vulnerable urbanizing parts of East-Africa.  _Current Opinion in Environmental Science & Health_, 5: 73–78.
  * Kratzert et al. (2019) Kratzert, F.; Herrnegger, M.; Klotz, D.; Hochreiter, S.; and Klambauer, G. 2019\.  NeuralHydrology–interpreting LSTMs in hydrology.  In _Explainable AI: Interpreting, explaining and visualizing deep learning_ , 347–362. Springer.
  * Little and Rubin (2019) Little, R. J.; and Rubin, D. B. 2019.  _Statistical analysis with missing data_ , volume 793.  John Wiley & Sons.
  * Maraun (2013) Maraun, D. 2013.  Bias correction, quantile mapping, and downscaling: Revisiting the inflation issue.  _Journal of Climate_ , 26(6): 2137–2143.
  * Mwale, Adeloye, and Rustum (2012) Mwale, F. D.; Adeloye, A. J.; and Rustum, R. 2012.  Infilling of missing rainfall and streamflow data in the Shire River basin, Malawi–A self organizing map approach.  _Physics and Chemistry of the Earth, Parts A/B/C_ , 50: 34–43.
  * Nash and Sutcliffe (1970) Nash, J.; and Sutcliffe, J. 1970.  River flow forecasting through conceptual models part I — A discussion of principles.  _Journal of Hydrology_ , 10(3): 282–290.
  * Norazizi and Deni (2019) Norazizi, N. A. A.; and Deni, S. M. 2019.  Comparison of artificial neural network (ANN) and other imputation methods in estimating missing rainfall data at Kuantan station.  In _International Conference on Soft Computing in Data Science_ , 298–306. Springer.
  * Pantanowitz and Marwala (2009) Pantanowitz, A.; and Marwala, T. 2009.  Missing data imputation through the use of the random forest algorithm.  In _Advances in computational intelligence_ , 53–62. Springer.
  * Pastén-Zapata et al. (2020) Pastén-Zapata, E.; Jones, J. M.; Moggridge, H.; and Widmann, M. 2020.  Evaluation of the performance of Euro-CORDEX Regional Climate Models for assessing hydrological climate change impacts in Great Britain: A comparison of different spatial resolutions and quantile mapping bias correction methods.  _Journal of Hydrology_ , 584: 124653.
  * Piani et al. (2010) Piani, C.; Weedon, G.; Best, M.; Gomes, S.; Viterbo, P.; Hagemann, S.; and Haerter, J. 2010.  Statistical bias correction of global simulated daily precipitation and temperature for the application of hydrological models.  _Journal of Hydrology_ , 395(3): 199–215.
  * Pigott (2001) Pigott, T. D. 2001.  A review of methods for missing data.  _Educational research and evaluation_ , 7(4): 353–383.
  * Rasmussen and Williams (2006) Rasmussen, C. E.; and Williams, C. K. I. 2006.  _Gaussian Processes for Machine Learning_.  The MIT Press.
  * Ringard, Seyler, and Linguet (2017) Ringard, J.; Seyler, F.; and Linguet, L. 2017.  A quantile mapping bias correction method based on hydroclimatic classification of the Guiana shield.  _Sensors_ , 17(6): 1413.
  * Rischard, Pillai, and McKinnon (2018) Rischard, M.; Pillai, N.; and McKinnon, K. A. 2018.  Bias correction in daily maximum and minimum temperature measurements through Gaussian process modeling.  _arXiv preprint arXiv:1805.10214_.
  * Sanchez Lozano et al. (2021) Sanchez Lozano, J.; Romero Bustamante, G.; Hales, R. C.; Nelson, E. J.; Williams, G. P.; Ames, D. P.; and Jones, N. L. 2021.  A Streamflow Bias Correction and Performance Evaluation Web Application for GEOGloWS ECMWF Streamflow Services.  _Hydrology_ , 8(2): 71.
  * Snow (2015) Snow, A. D. 2015.  _A New Global Forecasting Model to Produce High-Resolution Stream Forecasts_.  Brigham Young University.
  * Souffront Alcantara et al. (2019) Souffront Alcantara, M. A.; Nelson, E. J.; Shakya, K.; Edwards, C.; Roberts, W.; Krewson, C.; Ames, D. P.; Jones, N. L.; and Gutierrez, A. 2019.  Hydrologic modeling as a service (HMaaS): a new approach to address hydroinformatic challenges in developing countries.  _Frontiers in Environmental Science_ , 7: 158.
  * Tencaliec et al. (2015) Tencaliec, P.; Favre, A.-C.; Prieur, C.; and Mathevet, T. 2015.  Reconstruction of missing daily streamflow data using dynamic regression models.  _Water Resources Research_ , 51(12): 9447–9463.
  * Van Buuren (2018) Van Buuren, S. 2018.  _Flexible imputation of missing data_.  CRC press.
  * Wang, Zhang, and Villarini (2021) Wang, C.; Zhang, W.; and Villarini, G. 2021.  On the use of convolutional Gaussian processes to improve the seasonal forecasting of precipitation and temperature.  _Journal of Hydrology_ , 593: 125862.
  * Xu and Liang (2021) Xu, T.; and Liang, F. 2021.  Machine learning for hydrologic sciences: An introductory overview.  _Wiley Interdisciplinary Reviews: Water_ , 8(5): e1533.