Potential of the Coupled WRF/WRF-Hydro Modeling System for Flood Forecasting in the Ouémé River (West Africa)

Abstract

Since the beginning of the 2000s, most of the West-African countries, particularly Benin, have experienced an increased frequency of extreme flood events. In this study, we focus on the case of the Ouémé river basin in Benin. To investigate flood events in this basin for early warning, the coupled atmosphere–hydrology model system WRF-Hydro is used, and analyzed for the period 2008–2010. Such a coupled model allows exploration of the contribution of atmospheric components into the flood event, and its ability to simulate and predict accurate streamflow. The potential of WRF-Hydro to correctly simulate streamflow in the Ouémé river basin is assessed by forcing the model with operational analysis datasets from the European Centre for Medium-Range Weather Forecasts (ECMWF). Atmospheric and land surface processes are resolved at a spatial resolution of 5 km. The additional surface and subsurface water flow routing are computed at a resolution of 500 m. Key parameters of the hydrological module of WRF-Hydro are calibrated offline and tested online with the coupled WRF-Hydro. The uncertainty of atmospheric modeling on coupled results is assessed with the stochastic kinetic energy backscatter scheme (SKEBS). WRF-Hydro is able to simulate the discharge in the Ouémé river in offline and fully coupled modes with a Kling–Gupta efficiency (KGE) around 0.70 and 0.76, respectively. In the fully coupled mode, the model captures the flood event that occurred in 2010. A stochastic perturbation ensemble of ten members for three rain seasons shows that the coupled model performance in terms of KGE ranges from 0.14 to 0.79. Additionally, an assessment of the soil moisture has been developed. This ability to realistically reproduce observed discharge in the Ouémé river basin demonstrates the potential of the coupled WRF-Hydro modeling system for future flood forecasting applications.

1. Introduction

In its fifth report, the Intergovernmental Panel on Climate Change (IPCC) stresses the increment of the number of extreme weather events for the 21st century due to climate change [1]. Tropical countries of West Africa are threatened by climatic hazards, such as droughts, floods, high winds, the elevation of the sea level, etc. Droughts and floods are the most important in terms of damages and impacts. Countries in West Africa such as Benin, Burkina Faso, Cote d’Ivoire, Niger, Senegal, and Togo have suffered from catastrophic floods [2] with severe consequences including loss of life, property, and damage of property. Floods happen in geomorphologically low-lying areas, such as rivers, basins, lakes, etc., and occur often with heavy rain associated with severe thunderstorms, hurricanes, tropical storms, or meltwater from ice or snow flowing over ice sheets or snowfields [3]. Ouémé river, the main watershed in the Benin Republic, is used for irrigation, hydroelectricity, transport, and water resource management. Its utility led to the construction of hydrological structures; however, this construction exposes the residents on the riverfronts to flooding risks. In 2010, Ouémé river has recorded an overflow driving the country (Benin) to experience dramatic flooding, which affected 680,000 people, leading to 43 deaths [4] over the length of the Ouémé river (at Bétérou, Savè, and Bonou outlets).. To prevent these damages during the rainy season, the ideal solution proposed by Benin is the improvement of national and decentralized capacities to better manage and plan for extreme weather events and climate-related risks through the development of an early warning system (EWS) and better climate information (IC). Therefore, citizens are earlier informed about the level of the risks based on color patterns, which depend on the water heights at the level of the water gauges. The flooding water heights limit is defined based on historical rainfall episode statistics. In addition, EWS warns when there is a risk that the rise of the water level upward of the river could spread downstream.

Such destructive floods make it necessary to improve hydrological forecasts to reduce the vulnerability of regional and local communities. Identifying the drivers of change in flood regimes in West African watersheds is a complex task due to the heterogeneity of the region and the changing hydrological functioning of watersheds as a result of human activity. Ref. [5] found that there is a significant trend both in terms of flooding magnitude and frequency in West Africa, with two main patterns: the Sahelian showed increasing flood trends, whilst some Sudanian areas presented decreasing flood trends. Since IPCC [1] announced and highlighted an increase in the frequency of extreme climate events, and also the potentially significant increase of their magnitudes due to climate change, identification of accurate flood forecasting methods is required. The introduction of flood forecasting and early warning systems, the building of population awareness and preparedness, urban planning, discouraging human settlements in flood-prone areas, and the development of local institutional capacities, as well as forecasting methods development, should be appropriate actions to plan.

To determine an accurate flood forecasting method, several studies applied hydrological modeling based on conceptual, semi-distributed, physical-based, data-driven or machine leaning methods. The conceptual methods adopted in some previous studies [2,6,7,8] have been proved to be useful for hydrological modeling options because most are successful for rainfall-runoff simulation and flood modelling. The data-driven approach, based on the empirical relationship amongst datasets, which is usually driven through statistical constraints [9,10,11], has been investigated; despite their advantages (no need of physical equations and parameters or catchment characteristics), they present some limitations such as static models which cannot evaluate the change (e.g., land-use/land-cover change). With the increasing development in the machine leaning field, numerous methodologies are applied in operational flooding forecasting [12,13]. Thus, some studies [14,15,16] have proven that the machine learning-based model, performs better when it is compared with the physical-based, conceptual, and data-driven models. These methodologies do not take into account the influence of both atmosphere and land surface complexity in the performance of the forecasting.

The operational global weather forecast centers routinely provide relatively coarse precipitation forecasts with resolutions of 16–27 km (e.g., ref. [17]). Because these kinds of forecasts could not provide necessary details of complex, intense precipitation structures led by the mesoscale orography, land-surface heterogeneities, and land-water contrasts, ref. [18,19] used the Weather Research and Forecast (WRF) model to provide high-resolution precipitation forecasts at a 1.3–4 km horizontal resolution. They found that the WRF model was able to provide precipitation forecast both in terms of amount and spatial distribution. Additionally, ref. [17] used the output from the WRF model (over the upper-Jordan River basin) as forcing to run the Hydrological Model for Karst Environment (HYMKE). Similar experiments were set in different areas [20,21,22,23], and the authors showed that the precipitation estimates with WRF are able to reproduce the spatial distribution of precipitation, but may underestimate the magnitude of the heavy precipitation events when compared to rain gauges or global observation datasets.

Refs. [24,25,26,27,28,29,30] have shown the benefits of choosing coupled atmospheric-land-surface models for estimating temperature and precipitation in different areas. Ref. [31] investigated and evaluated simulations with a fully coupled atmospheric–hydrological model (WRF-HMS) and an uncoupled model for various meteorological variables, and demonstrated the better performance of the fully coupled model against the uncoupled. Another study led by [32] compared a one-way forced implementation of the WRF-Hydro system to a fully two-way coupled instance of WRF/WRF-Hydro in order to assess the impact of two-way coupling on simulated precipitation and streamflow. They found that the two setups performed well for the precipitation, but the matching to observed data was higher for the two-way coupled WRF/WRF-Hydro simulation in terms of statistical performance criteria.

Assessing the advantages and limitations of one-way vs. two-way coupled modeling systems for flood prediction over the Ayalon basin (Israel), ref. [17] used both a hydrological model (Hydrological Engineering Center-Hydrological Modeling System, HEC-HMS) and the WRF-Hydro modeling system. The models were forced by observed, interpolated precipitation from rain-gauges within the Ayalon basin, and with modeled precipitation from the WRF atmospheric model. Comparing simulations with the one-way coupled WRF model and the two-way coupled WRF/WRF-Hydro modeling system, they found that the use of two-way atmospheric-hydrological coupling has the potential to enhance precipitation simulation and, therefore, hydrological forecasts for early flood warning applications. It is important to acknowledge here that the model uncertainty in simulating precipitation, as well as streamflow, can be relatively high. Model uncertainty can, for example, be evaluated with a model ensemble (e.g., [33,34]). Ref. [35] developed the stochastic kinetic energy backscatter scheme (SKEBS) within the WRF model in order to generate such an ensemble.

The present study attempts to develop a model system that is able to simulate flood events in the Ouémé river basin (Benin, West Africa). According to successful examples in the literature, we use the fully coupled atmospheric-hydrological modeling system WRF/WRF-Hydro for this purpose. Section 2 gives the characteristics of the study area, different sources of the dataset used, and the modeling approaches. Results and discussion are provided in Section 3, whilst Section 4 is dedicated to a summary and conclusion.

2. Materials and Methods

2.1. Study Area and Observation Datasets

The study area is in West Africa, located between latitudes 0° N and 18° N and longitude 7° W and 12° E, as displayed in Figure 1. This also constitutes the setup of the WRF domain. The region is bordered in the south by the Gulf of Guinea, and in the north by Mali and Niger. Nigeria’s highlands form the eastern boundary, while Mauritania, Mali, and Ghana form the western limit. The annual mean temperature is about 18 °C, but the monthly mean can be more than 30 °C over the southern part of the Sahara. The rainfall pattern over this region is mostly affected by ocean currents and local features such as topography. In terms of climatic zones, West Africa is characterized by mainly three different regions: the first region covers the Sahel and is characterized as a semi-arid zone located from western Senegal to eastern Sudan between 12° N and 20° N. The second zone is the Sudano-Sahelian region, and the third zone comprises the Guinea coast, which is characterized by a bimodal mode driven by the inter-tropical divergence (ITD). The basin is located in the Benin Republic.

The Republic of Benin is in the inter-tropical zone (between 06°10′ N and 12°25′ N), which has a wet and dry tropical climate [2]. It contains the rivers (Savè and Bétérou) on which we focus in this study (see Figure 1). The Ouémé catchment at the Savè (resp. Bétérou: inner-catchment to Savè) outlet covers an area of 24.800 km² (resp. 10.475 km²). It is located between 7°58′–10°12′ N and 1°35′–3°05′ E, and represents 47% of the whole Ouémé river (ref. [36]). The seasons correspond to the periods of dominance of the wet tropical continental air masses. The seasonal distribution of rainfall follows the direction of the ITD and varies almost proportionally with distance from the coast. Therefore, Bétérou has a unimodal precipitation regime (May to October), whilst the southern part of the Savè catchment has a transitional regime (April and October, with a short dry period in August). The average annual rainfall between 1960–2007 is 1200 mm at the Bétérou rainfall station, and 1100 mm at Savè. The flow dynamic is characterized by a high discharge during the rainy season. The maximum flow between May and September over the period 1960–2007 is in the order of 270 m³/s at Bétérou and 480 m³/s at the Savè outlet. From November to May, almost all the rivers dry up with averages of low flows of 5 m³/s at Savè and 2 m³/s at Bétérou. The annual mean temperature range is between 24 °C and 33 °C.

Our assessment is focused on the year 2010, known as the year when the Ouémé river experienced dramatic flooding as a case study. Ref. [2] showed that the maximum values of discharge recorded during the period 1989–2009 is less than 1400 m³/s at Savè, and 650 m³/s at Bétérou. Analysis of station data for the period 1960–2007 showed that the peaks of discharge at Savè (resp. Bétérou) are about 910 (resp. 470), 1067 (resp. 560), and 1200 (resp. 640) m³/s, respectively, for a 5-, 10-, and 20-year return period. The discharge and precipitation station data used in this study were collected over Savè and Bétérou outlets of the Ouémé river basin. The three-hourly satellite estimates of the Tropical Rainfall Measuring Mission (TRMM, 3B42 v7 derived daily at 0.25° horizontal resolution, 1998–near-present; ref. [37]) dataset, and the daily Climate Hazards Group Infrared Precipitation with Stations (CHIRPS; chirps v2.0 at 0.05° horizontal resolution; 1981–near-present; ref. [38]) is used for model evaluation.

2.2. Experimental Design of the Weather Research and Forecasting (WRF) and Coupled WRF/WRF-Hydro Model

The advanced Weather Research and Forecasting (WRF, [39,40]), developed by the U.S. National Center for Atmospheric Research (NCAR) version 3.7.1, is used for both the WRF-only and fully coupled WRF/WRF-Hydro modeling in the study. The fully coupled WRF/WRF-Hydro, hereafter referred to as WRF-H, is a non-hydrostatic, mesoscale numerical weather prediction and atmospheric simulation system.

Table 1. Experimental details of the atmosphere model, WRF-only and WRF-H.

Subject	Option	Reference
Driving data	European Centre for Medium-Range Weather Forecasts (ECMWF) operational analysis	ECMWF
Horizontal resolution	5 km
Horizontal grid	400 × 400
Integration time step	30 s
Projection resolution	Mercator
Vertical discretization	50 layers
Output interval	24 h for WRF, 30 days for WRF/WRF-Hydro
Simulation period	1st January 2008–31st December 2010
Pression top	10 hPa
Microphysics scheme	Single Moment Microphysics class 5 (WSM5)	[43]
Longwave radiation	Rapid Radiative Transfer Model (RRTM)	[44]
Shortwave radiation	Dudhia	[45]
Planetary boundary lager	Asymmetric Convection Model (ACM2)	[46,47]
Land use	Moderate Resolution Imaging Spectroradiometer (MODIS)	[42]
Land surface scheme	Noah Land Surface Model (LSM)	[25]

summaries the different physics schemes and experimental details. The model setup uses one domain at 5 km spatial resolution, covering the area 7° W–12° E, 0°–18° N, and 400 × 400 grid points, with 30 s as numerical simulation time step. The vertical structure of the domain consists of 50 levels, up to a 10 hPa pressure top. The WRF-Hydro model is an extension package of a hydrological component to WRF. This supplementary model allows the lateral representation of the hydrological condition at the land surface [41]. The option of land use categories “Moderate Resolution Imaging Spectroradiometer (MODIS, 20 classes; ref. [42])” is selected. The Noah LSM model ref. [25] is used as the column land surface physics model.

For purposes of hydrometeorological simulations with WRF-H, the WRF domain is additionally coupled with routing processes at 500 m resolution for 4000 × 4000 grid points in the east-west and north-south directions. The fully coupled mode simulations are performed for three years, from January 2008 to December 2010, with January–February 2008 as the spin-up period. The driving data is the operational analysis dataset from the European Centre for Medium-Range Weather Forecasts (ECMWF), which provides the initial and lateral boundary conditions at 0.125° horizontal resolution and 6 h time intervals. Both WRF-H and WRF-only components of the coupled modeling system share the same physics parameterizations in Table 1.

The WRF is the atmospheric model which solves the equations of atmospheric motions. The atmospheric state variables are the wind components, pressure, temperature, humidity and hydrometeor mixing ratios, etc. The lower boundary of the atmospheric variables is forced with a land surface model. In the case of WRF-Hydro, a complex land surface model, with an explicit description of river streamflow, is employed. The uncoupled WRF-Hydro model consists of a variety of parameters (e.g., ref. [48]), which usually require calibration. Since the aim of the research is to evaluate the performance of WRF-H to simulate discharge, and therefore analyze its predicting skills about floods, the calibration is performed based on discharge at the Savè catchment outlet (Figure 1).

2.3. Calibration of WRF-Hydro in Offline Mode

For calibrating the model WRF-Hydro 4.0, we focus on selected sensitive parameters highlighted in previous studies, such as REFKDT, SLOPE, RETDEPRTFAC, OVROUGHRTFAC, and MannN. Applying a stepwise approach, following previous WRF-Hydro studies [17,32,49,50], we first focus on the parameters controlling the total water volume, namely infiltration factor (REFKDT) and surface retention depth (RETDEPRT, see

Table 2. Validation of the performance of simulated WRF-Hydro and observed discharge at Savè, with parameters such as the infiltration-runoff parameter REFKDT, the retention factor RETDEPRTFAC, the SLOPE, the overland flow roughness scaling factor OVROUGHTFAC, and the Manning’s roughness coefficients, MannN. Experiments in bold show the optimum parameters based on the objective criteria (Corr, NSE, and KGE).

REFKDT
Range	0.1	0.5	0.8	1.0	1.5	2.0	3.0	3.5	4.5
NSE	0.29	0.33	0.38	0.46	0.52	0.34	0.20	0.17	0.11
Corr	0.64	0.67	0.63	0.63	0.58	0.60	0.61	0.64	0.60
KGE	0.38	0.39	0.41	0.47	0.49	0.32	0.09	0.07	0.07
RETDEPRTFAC
Range	0.0	1.0	2.0	3.0	4.0	5.0	6.0	7.0	8.0
NSE	0.50	0.52	0.49	0.49	0.49	0.49	0.49	0.49	0.49
Corr	0.58	0.58	0.59	0.59	0.59	0.59	0.59	0.58	0.58
KGE	0.48	0.49	0.49	0.49	0.48	0.46	0.47	0.46	0.46
SLOPE
Range	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9
NSE	0.52	0.61	0.60	0.50	0.37	0.20	0.12	−0.03	−0.23
Corr	0.58	0.62	0.65	0.65	0.66	0.61	0.63	0.55	0.61
KGE	0.49	0.56	0.52	0.39	0.33	0.24	0.19	0.11	0.04
OVROUGHRTFAC
Range	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1.0
NSE	0.60	0.60	0.60	0.61	0.63	0.62	0.62	0.61	0.61
Corr	0.67	0.65	0.61	0.60	0.66	0.66	0.64	0.64	0.62
KGE	0.56	0.56	0.57	0.60	0.60	0.60	0.59	0.56	0.56
MannN
Range	0.2	0.4	0.6	0.8	1.0	1.4	1.6	1.8	2.0
NSE	0.42	0.44	0.50	0.56	0.61	0.62	0.66	0.62	0.60
Corr	0.66	0.63	0.73	0.68	0.70	0.72	0.67	0.69	0.67
KGE	0.36	0.37	0.41	0.48	0.50	0.57	0.63	0.60	0.61

). It is noted that REFKDT is a tunable parameter that significantly impacts surface infiltration, and hence the partitioning of total runoff into the surface and subsurface runoff; increasing REFKDT decreases surface runoff. Since there is not a historical range to estimate these parameters over the interested domain, the study calculates them from 0.1 to 10, with 0.1 increments. The second step of the calibration is to evaluate the coefficient governing deep drainage (SLOPE); the same method used in the case of REFKDT and RETDEPRTFAC for selecting the optimum value is applied, by testing values from 0.1 to 1.0, with 0.1 as the increment. The adjustment of the roughness parameter, which controls the overland flow, is performed from the default value to the optimum one.

In addition, sensitivity tests are done on the surface (OVROUGHRTFAC) and channel roughness parameter (MannN), which controls the shape of the hydrograph. The efficiency criteria—Kling–Gupta efficiency (KGE), Correlation coefficient (Corr), Nash- Sutcliffe Efficiency (NSE) and BIAS are used to evaluate the model performance within the calibration process.

$$Corr=\frac{n\left(\sum xy\right)-\left(\sum x\right)\left(\sum y\right)}{\sqrt{\left[n\sum x^2-{\left(\sum x\right)}^2\right]\left[n\sum y^2-{\left(\sum y\right)}^2\right]}}$$

(1)

$$NSE=1-\left[\frac{{\sum }_{i=1}^n{\left(Y_i{}^{obs}-Y_i{}^{sim}\right)}^2}{{\sum }_{i=1}^n{\left(Y_i{}^{obs}-Y_i{}^{mean}\right)}^2}\right]$$

(2)

$$KGE=1-\sqrt{{\left(r-1\right)}^2+{\left(\frac{{\sigma }_{sim}}{{\sigma }_{obs}}-1\right)}^2+{\left(\frac{{\mu }_{sim}}{{\mu }_{obs}}-1\right)}^2}$$

(3)

$$BIAS=\frac{{\sum }_{i=1}^n\left(Y_i{}^{sim}-Y_i{}^{obs}\right)}{{\sum }_{i=1}^n{Y}_i{}^{obs}}\times 100$$

(4)

where $x$ is the observations and $y$ is the simulations; $Y_i{}^{obs}$ is the ith observation for the constituent being evaluated, $Y_i{}^{sim}$ is the ith simulated value for the constituent being evaluated, $Y_i{}^{mean}$ is the mean of observed data for the constituent being evaluated, and $n$ is the total number of observations [51]. For the KGE, $r$ is the linear correlation between observations and simulations, ${\sigma }_{obs}$ is the standard deviation in observations, ${\sigma }_{sim}$ is the standard deviation in simulations, ${\mu }_{sim}$ is the simulation mean, and ${\mu }_{obs}$ is the observation mean.

In order to harmonize the uncoupled and coupled setups, the uncoupled simulations use the same time step as the WRF-only and WRF-H simulations (30 s). The calibration of the model is performed using hourly dataset input, and the focus is on the performance skill in reproducing daily discharge in the sub-catchments. One-year calibration is considered sufficient to evaluate the basic parameter sensitivities [32]. WRF-Hydro is calibrated on period P1 (2008) and validated on period P2 (2009–2010), where P1 and P2 are the shared periods containing the whole study period, named P (2008–2010).

2.4. Evaluation of Model Uncertainty with the Stochastic Kinetic Energy Backscatter Scheme

The stochastic kinetic energy backscatter scheme (SKEBS; [35,52,53]), which primarily acts on the dynamical tendencies at the lateral boundaries, is activated into WRF-H for the fully coupled simulation (WRF-H-SKEBS). The SKEBS technique provides several advantages over perturbation techniques that only perturb the initial state. The method aims to represent model uncertainties associated with scale interactions that take place in the real atmosphere but are absent in a truncated numerical model [54]. SKEBS perturbs the model fields by adding random and amplitude perturbations (noise) to the horizontal wind and potential temperature tendency equations at the lateral boundaries for each time step [55]. An ensemble of 10 members using the WRF-H-SKEBS model is generated for each rain season of the period P (2008–2010).

3. Results and Discussion

3.1. Evaluation of WRF-Only Precipitation

The WRF model (no coupled, referee as WRF-only) is run over the domain displayed in Figure 1 in order to generate atmospheric input data for the calibration of the supplement model WRF-Hydro. Figure 2 shows the comparison between the simulated WRF-only and observed dataset. The weekly WRF precipitation for the Savè catchment is relatively close to that derived from CHIRPS (Figure 2a) and TRMM (Figure 2b), with the mean coefficient of determination (R²) equal to 0.64, and 0.59, respectively. The agreement between the two observed datasets (CHIRPS and TRMM) is about 0.87 for R².

Since the investigation is on the potential of WRF-H (model WRF-Hydro coupled with WRF) for flood prediction, the available WRF-only precipitation at the highest spatiotemporal resolution is used to force the uncoupled WRF-Hydro model. In particular, the hourly output of WRF at 5 km spatial resolution is used as meteorological forcing data, containing necessary variables such as incoming shortwave radiation (W/m²), incoming longwave radiation (W/m²), specific humidity (kg/kg), air temperature (K), surface pressure (Pa), u and v components of near-surface wind (m/s), and liquid water precipitation rate (mm/s). The meteorological forcing data needed by the Noah LSM (land surface hydrological modeling system) are prepared as hourly gridded data. The Noah LSM static data (topography, land cover, soil type) are too coarse for a WRF-Hydro application. Therefore, datasets from the Shuttle Elevation Derivatives at Multiple Scales (HydroSHEDS) database [56], high-resolution topography and channel network) are considered to accurately route water across the landscape through overland, subsurface, or channel flow.

3.2. Calibration and Evaluation of WRF-Hydro Offline

The Noah LSM coupled with WRF-Hydro evaluates surface water depth change h(m) as a rate over time of infiltration excess [25]:

$$\frac{\partial h}{\partial t}=\frac{\partial P_d}{\partial t}\left\{1-\frac{\left[{\sum }_{i=1}^4\Delta Z_i\left({\theta }_s-{\theta }_i\right)\right]\left[1-\exp \left(-k\frac{K_s}{K_{ref}}\frac{{\delta }_t}{86400}\right)\right]}{P_d+\left[{\sum }_{i=1}^4\Delta Z_i\left({\theta }_s-{\theta }_i\right)\right]\left[1-\exp \left(-k\frac{K_s}{K_{ref}}\frac{{\delta }_t}{86400}\right)\right]}\right\}$$

(5)

where $P_d$(m) is the precipitation that is not intercepted by the canopy; $\Delta Z_i\left(m\right)$ is the soil layer $i$ in Noah LSM, ${\theta }_i$ is the volumetric water content (soil moisture) of the soil layer $i$, ${\theta }_s$ is the saturated soil moisture, and the saturated hydraulic conductivity is $K_s \left(m/s\right)$; both ${\theta }_s$ and $K_s$ varied in terms of the soil texture. The reference saturated hydraulic conductivity $K_{ref}$ is equal to $2\times {10}^{-6}\left(m/s\right)$.

The calibration of the model (WRF-Hydro) is performed in offline mode for the Savè catchment with parameters REFKDT, SLOPE, RETDEPRTFAC, OVROUGHRTFAC, and MannN, using both KGE and Corr as efficiency criteria. Efficiency criteria results of the optimization of those parameters are listed in Table 2. The first parameter evaluated is one of the parameters controlling the total water volume (the REFKDT), whose values range from 0.1 to 10.0, with 0.1 increments. The default simulation (REFKDT = 3.0) realized had shown an underestimation of the observed streamflow because the REFKDT controls the infiltration of the surface water, and its value shall be reduced to disable many infiltrations. As illustrated with [49] for the case of the Sissili in West Africa, and [48] in Kenya (East Africa), we found that the model discharge performance is highly sensitive to parameter REFKDT. For the case in study, the REFKDT = 1.5 performs better than using the default value (3.0) with statistics (NSE = 0.52, KGE = 0.49, and Corr = 0.58). The second parameter assessed is the RETDEPRTFAC, which adjusts the initial retention depth in the model. Scaling factors (RETDEPRTFAC) between 0 and 10 with 1.0 increments show that the modeled discharge remains the same regardless of the change in RETDEPRTFAC, as the scores are very close to each other.

As mentioned by [41], increases in the RETDEPRTFAC on channel pixels can encourage more local infiltration near the river channel leading to wetter soils. Therefore, the default RETDEPRTFAC = 1.0 seems to provide better results and will be considered such for the next calibration steps. In the case of the present study, the uncoupled model WRF-Hydro is also very sensitive to the parameter SLOPE. The optimized value of the parameter SLOPE is 0.2. The parameters controlling the hydrograph shape are also investigated. As illustrated by [50], the surface roughness (OVROUGHRT) plays an important role in transmitting infiltration excess water to channel networks and is calibrated in WRF-Hydro using a scaling factor (OVROUGHRTFAC) between 0.2 and 1.0, with 0.1 increments. Considering the correlation coefficient (Corr), the NSE and the KGE statistics, and the hydrograph shape match between simulated and observed hydrographs at Savè, the scaling factor value of OVROUGHRTFAC = 0.6 is judged as the best to fit the simulated hydrograph to the observed hydrograph. Furthermore, the calibrated Manning’s coefficient (MannN) for the river channel routing used are set as 1.75 for stream order 1, 1.70 for stream order 2, 1.65 for stream order 3, 1.60 for stream order 4, and 1.55 for stream order 5. It can be seen that the observed discharge hydrograph at Savè is reasonably well reproduced, with KGE = 0.63, and Corr = 0.67, for the calibration period March-December 2008 (Figure 3a).

The calibrated model is evaluated offline for the period P2: 2009–2010 (Figure 3b). The above-mentioned efficiency criteria allow us to evaluate the performance of the model. It can be noticed that it simulates the trend and peaks of the observed discharge relatively well, even slightly better in comparison to the calibration period, with model efficiencies KGE of 0.86 and Corr of 0.87. This enhanced performance for the validation period is related to the much higher discharge peak in 2010, i.e., the flooding year, which is relatively well reproduced by the model evaluated offline. Globally, for the simulation period P (2008–2010), WRH-Hydro in offline mode is able to simulate discharge, with KGE and Corr equal to 0.70 and 0.74.

The robustness of the calibrated model in the uncoupled module is also evaluated over the Bétérou catchment (an inner-domain of the Savè basin). Figure 4 shows the best-fit performance of the model in simulating the shape of the observed hydrograph at Bétérou. It reproduces the peaks of discharge relatively well, and the trend of the shape of the observed hydrograph is also well reproduced. The performances of the model are revealed by the statistics KGE = 0.74, Corr = 0.85, and NSE = 0.68 for the calibration period over the Bétérou catchment for the uncoupled module simulation.

3.3. Evaluation of the Coupled Model WRF-H

The calibrated (optimized) model parameters are used for the online WRF-Hydro model (referred to as WRF-H), to assess the performance of the calibrated WRF-H to simulate discharge and precipitation in the study area.

3.3.1. Precipitation Simulations

The agreement skills of WRF-H are evaluated for precipitation both temporally and spatially for the rainy season of the period 2008–2010. Figure 5 exposes a comparison between the weekly precipitations from WRF-only, WRF-H, and observed datasets. The coefficient of determination between the WRF-H and WRF-only, is equal to 0.88 (Figure 5a). This shows clearly that WRF-only and WRF-H simulate precipitation differently, which has been already illustrated by [17,32,57,58]. Figure 5b,c compare WRF-H with CHIRPS and TRMM, and it illustrates a good agreement between these datasets. The slightly better agreement of CHIRPS (compared to TRMM) with WRF-only (Figure 2b) and WRF-H (Figure 5b) could be explained by the high resolution of both CHIRPS and WRF-H precipitation. As seen previously (Figure 2d), CHIRPS and TRMM precipitations’ have the same trend, only differ in their magnitudes; the WRF-only depicts well the observations for the whole season. Indeed, the comparison of Figure 2 and Figure 5 illustrates that WRF-H performs slightly better than WRF-only in terms of weekly precipitations. Overall, Figure 5d shows that WRF-H follows well in terms of the trend, the WRF-only. Additionally, Figure 5d enhances these results with Corr = 0.68 between WRF-H and CHIRPS, against Corr = 0.59 between WRF-H and TRMM. Despite the good performance of models, they present some incertitude wherever they are applied. Ref. [59] showed that the high variability of precipitation in West Africa results from a large uncertainty in WRF simulations. This uncertainty is investigated in detail in the case of WRF-H in Section 3.4 by modifying boundary conditions with a stochastic perturbation.

The monthly trend of precipitation during the West African Monsoon (WAM) period (June–September: JJAS) has been investigated at the Ouémé–Savè outlet (Savè). It shows that the monthly cumulative trend follows observations shape but varies in magnitude, as illustrated in Figure 2 and Figure 5. Figure 6 shows the maximum monthly cumulative rainfall at Savè in August, and while comparing the three patterns, it shows that the monthly cumulative precipitation records during the two last months (August and September) in 2010 are higher than both 2008 and 2009, expressing the heavy rainfall recorded in this particular year. To observe more detail of the distribution of precipitation over the focusing domain for either WRF-only and WRF-H, as well as CHIRPS, Figure 7 presents the spatial distribution of precipitation during the WAM (JJAS) of the flooding year (2010).

Figure 7a,b depict the WRF-only simulations and WRF-H precipitations, respectively, while the difference between the two models in Figure 7c shows that WRF-H either underestimated or overestimated simulation in comparison to WRF-only, depending on the location (as illustrated [31]). Figure 7d depicts the spatial distribution of CHIRPS over the focusing domain; its color bar is as in Figure 7a,b. The mean precipitation for the domain is about 864 mm for WRF-only (Figure 7a), and 947 mm for WRF-H (Figure 7b). This means that WRF-H increases the simulated precipitation from WRF-only by about 1%. The observed precipitation in this domain is about 817 mm (i.e., less than simulated precipitation from both WRF-only and WRF-H). Similar results are obtained for the Savè catchment, with seasonal spatial-averaged precipitation of about 1049 mm for WRF-H, 998 mm for WRF-only, and 977 mm for CHIRPS.

As the model performed globally well precipitations both for Savè and Bétérou at temporal scale and spatial distribution, the following analysis focused on those two basins of the Oueme river by investigating the temporal cumulative variability of precipitation. Therefore, Figure 8 and Figure S1 present, respectively, the cumulative precipitations derived from WRF-H, and the satellite datasets CHIRPS and TRMM at Beterou and Savè. They attempt to assess WRF-H behavior about the cumulative precipitation over the study period (2008–2010). A very good reproductivity of the observed dataset by WRF-H is noticed, and also the clear rainy season period is well captured by the model. It can also be noticed that the total seasonal precipitation is based on the rainfall from April to October over the basins. The cumulative seasonal amount (between April and October) is indicative of the average of annual precipitation for a given year in the region.

According to Figure 8 and Figure S1, we can make the assumption that there is no or negligible rainfall recorded in the dry season. At Bétérou (Figure 8), the model underestimates the cumulative precipitation compared to TRMM and CHIRPS, albeit showing a good capture of seasonal variability over the years. The respective total three-yearly amounts are WRF-H = 3921 mm, TRMM = 4077 mm, and CHIRPS = 4000 mm, which confirms the consistency with the results discussed above that the WRF-H precipitation simulations are closer to the observed CHIRPS than TRMM. The bias of the cumulative precipitation is between −8.82% and 0.00% in regard to CHIRPS, and between −8.16% and 2.13% for TRMM. At Savè (with WRF-H = 3883 mm), conversely to analysis at Bétérou, the cumulative precipitation of the model fits well with TRMM (that value is about 3985 mm) rather than CHIRPS (value equal to 4020 mm). The error related to the evaluation of the cumulative precipitation is between −5.82% and 2.18% for TRMM, and about −7.00% and 0.00% for CHIRPS.

3.3.2. Discharge Simulations

Discharge results at Savè are displayed in Figure 9a, showing the daily time series of simulated (green) and observed (red) stream discharges, and related WRF-H precipitation (blue) for the period 2008–2010.

A good agreement can be seen between the observed and the simulated hydrographs, and an approximate good representation of the peaks of discharge, as well as hydrograph shapes, as quantified by the performance measures KGE and Corr, equal to 0.76 and 0.84, respectively. This better performance, in comparison to the offline simulation, could be explained by the time step of the meteorological data in the fully coupled mode, which is 30 s and not hourly as in offline mode. Since one objective of the study is to evaluate the performance of WRF-H to simulate the discharge, and therefore to predict potential floods, here we focus on the ability of the model in reproducing only the rainy seasons. We obtain for Savè, from Figure 9, a KGE equal to 0.22, 0.64, and 0.80 for the rainy seasons of 2008, 2009, and 2010, respectively, which gives solid information about the model’s simulation skills. It is noted that the model has a better performance in 2010. The robustness of the calibrated WRF-H over Savè is evaluated in a second catchment, i.e., in Bétérou (Savè’s inner-catchment).

Figure 9b shows that WRF-H here reproduces well the discharge trend, as well as the peaks, so that WRF-H can also be used successfully for this inner-catchment. Table S1 illustrates the discharges peaks obtain for the two basins during the three years. An evaluation from a recent work of ref. [2] over the same study domain revealed that the averages of annual rainfall between 1960–2007 are 1205 mm at the Bétérou rainfall station and 1098 mm at Savè. The dynamics of the flow are characterized by a high discharge during the rainy season. The maximum flow between May and September over the period 1960–2007 is in the order of 240 m³/s to 740 m³/s at Bétérou, and 1000 m³/s to 1750 m³/s at the Savè outlet. From November to May, almost all the rivers dry up, and the averages of low flows are about 5 m³/s at Savè and 2 m³/s at Bétérou.

WRF-H is able to capture the flood event which occurred in September–October 2010 over Savè, as well as over Bétérou. In particular, although the predicted highest discharge peak occurs earlier than in the observation at Savè and Bétérou. The second “weak” peak in 2010, which could amplify damage intensities of the flood in the study area, is also well reproduced. According to results from Figure 9, this second peak noticed for the discharge, should result from the highest precipitation simulated and observed in September 2010. The first important peak at Savè in 2010 is also reflected from the highest simulated precipitation of August 2010.

Figure 10 and Figure 11, respectively, over Bétérou and Savè show, on the one hand, the simulated and observed discharge at daily time step for the whole WRF-H evaluation period (2008–2010: Figure 10a,b; Figure 11a,b), and the flooding period (2010: Figure 10c; Figure 11c) based on linear regression with values of line 1:1 and coefficient of determination (R²); and on the other hand, the cumulative discharge of simulated and observed datasets. Figure 11a indicates that, at Bétérou, the WRF-H underestimates the low flows, whilst in contrast, at Savè (Figure 11a), it extracts the low flow relatively well. This underestimation of the low flow at Bétérou and Savè is related to the fact that in early August at the beginning of the rainy season over the sub-basins, the model started producing streamflow (2008 and 2009 on Figure 10a,b). Figure 9a justifies that the underestimation of the low flow at Savè (Figure 11) is due to the lags noticed between the simulated and observed “peaks”.

In general, WRF-H captures well the variability of the seasonal discharge for Savè and Bétérou, but overestimates the cumulative total discharge (Figure 10c and Figure 11c). These figures also indicate that the discharge is recorded from July to November, with a high increase in August–September. The total volume of discharge observed at Bétérou during the period 2008–2010 is 86,045 m³ against 100,182 m³ for the WRF-H simulation, with a rise of 16% as a difference. At Savè, the recorded total volume of water is about 186,475 m³, while the corresponding simulated discharge is approximately 220,367 m³, overestimating the observed discharge up to 18%.

3.4. Evaluation of the Soil Water Content

Figure 12 and Figure S2 show the daily averaged soil moisture extracted for soil depth between 0.5–2 cm [60,61] from the Climate Change Initiative (CCI) of the European Space Agency (ESA: https://www.esa.int/ESA (12 May 2018).; thereafter referred to as ${\theta }_{CCI}$), respectively, for Bétérou and Savè, and the corresponding volumetric soil water content in the first Noah soil layer between 0 and 10 cm (referred to as ${\theta }_{WRF-H}$). The ${\theta }_{WRF-H}$ is relatively high at the beginning of the simulation (January 2008: about 0.3 m³·m⁻³), but experiences a decreasing trend from the first months of the simulation, and reached approximatively the same value in February as those of other years (February 2009 and February 2010). This suggests that there is an excess of ${\theta }_{WRF-H}$ both at Bétérou and Savè at the initial time of the simulation, and that a two-month spin-up period appears to be sufficient for soil moisture in the first Noah LSM soil layer [33]. The ${\theta }_{CCI}$ values are globally lower than WRF-H simulations either for Bétérou or Savè during the rainy season and simulate reasonably well at the onset end the end of the season. However, despite the important bias, WRF-H correlated well with the observed.

Furthermore, Figure 13 and Figure S3 present the daily soil moisture in the four Noah LSM soil layers of the WRF-H simulation, respectively, for the basins Bétérou and Savè. As observed for the ${\theta }_{1WRH-H}$ (with soil layer between 0–10 cm) where only a two-month spin-up is enough to simulate the soil moisture, the difference between the two years is large for the soil layer 10–200 cm. For instance, for the soil moisture value of the second (10–40 cm) and the third (40–100 cm) soil layer, the stability is reached after the simulation of more than one year. Therefore, for this study case, we need more time for spin-up (more than a year) to adequately simulate the soil moisture between 10–200 cm.

3.5. Evaluation of Uncertainty of WRF-H

In order to evaluate the forecasting uncertainties of WRF-H (analyzed in Section 3.3), a stochastic kinetic energy backscatter scheme (SKEBS: [21,34,36]) is used and activated in WRF-H; it is referred to as WRF-Hydro-SKEBS. The purpose here is that the SKEBS approach adds random perturbations with prescribed spatial and temporal decorrelations. In particular, SKEBS produces perturbation into the lateral boundary conditions. The amplitude of the stochastic perturbations is chosen as the default in the WRF model. Therefore, we generated an ensemble of ten members for the evaluation of the uncertainties. Both stochastic physics and initial condition perturbations into WRF-Hydro-SKEBS result in an ensemble spread for the three rainy seasons over both Savè and Bétérou. Figure 14 shows that WRF-Hydro-SKEBS has a relatively large impact on precipitation and discharge results in the study region. The WRF-Hydro-SKEBS ensembles depicted well the discharge shape, but showed either an overprediction or underestimation of the discharge intensity. The assessment of the forecast uncertainty showed that the high variation of discharge amounts is associated with a speedy variation of the forecasted precipitation ensemble as demonstrated [62]. The ensemble also results in a large range of simulated discharge performance, as can be seen in Table S2. As in the case of the current study where the forecasting is performed under the convective-scale (5 km horizontal resolution), [63,64] concluded that SKEBS (at convective-scale) provides a balance between increased ensemble spread and forecast bias change. This demonstrates the sensitivity of WRF-H to lateral boundary perturbations and confirms the uncertainty of the model regarding discharge and precipitation simulations, which is of the uttermost importance for flood forecasting.

4. Summary and Conclusions

The present study explored the abilities of the fully coupled WRF-Hydro modeling system to simulate discharge and precipitation in the Ouémé river in West Africa. The main objective of the study was to determine an accurate flood forecasting tool that exhibits the interaction between atmosphere and land surface models, since the feedback of these interactions impact the precipitation dynamic [62,65]. Therefore, the supplementary hydrological extension of WRF model has been calibrated in offline mode for one year, and tested for two years using hourly outputs from WRF simulations. Optimized parameters from the calibration were used to perform the fully coupled WRF-Hydro model, which was used to investigate the forecasting performance skills of WRF over the Ouémé. We have found that the coupled modeling system WRF-H performed slightly better regarding the precipitation episodes than the WRF-only model (as already illustrated by [66] over the interested period (2008–2010)).

The evaluation of simulated precipitation showed its good performance skills, and provides confirmation about the uncertainty of WRF-H to simulate precipitation [43,46]. WRF-H also showed a good performance to simulate discharge, with a KGE equal to 0.76 for the period 2008–2010. The robustness of WRF-H has been assessed at Bétérou, an inner-catchment of the Ouémé river at Savè, where it provided a good agreement with respect to observed discharge, with a KGE equal to 0.66. Additionally, WRF-H was able to capture the flood event which occurred in 2010 over both Savè and Bétérou. Indeed, in the WRF-Hydro simulation in the fully coupled mode, the atmospheric and hydrological processes are simulated in a consistent way, which enhances the confidence in the results.

The uncertainty of predictability skills of WRF-H with respect to discharge in the Ouémé river at Savè was performed with an ensemble of ten members using a random perturbation scheme. Results demonstrated the large sensitivity of simulated discharge to perturbations introduced into the atmosphere. In summary, WRF-H is considered as a suitable model for evaluating discharge prediction uncertainties in the Ouémé river, and we encourage the implementation of the model for further basins in West Africa. Further studies are needed to assess how such WRF-H predicted discharges uncertainties could be used to improve an EWS in West-African regions.

Fully coupled WRF-H allows simulation of streamflow from global meteorological data. As such, coupled WRF-H can easily be applied to any region of the world, which is an advantage of the method used. Nevertheless, the application of WRF-H to another region may require a different calibration to optimize the quality of the modeled streamflow results over the concerned area. Additionally, the global forcing dataset to drive the model must be available and provided. Moreover, the model has some shortcomings; the human impacts (factors) on streamflow, such as irrigation and dam managements, are not yet considered in WRF-Hydro, which we see as a future and necessary model improvement for flood risk management purposes.