Figure 2 . Study area: Purus River Basin. (a) drainage network
(in blue), location of the discharge gauge (Canutama, triangle in red),
tracks of the spatial altimetry mission Jason 2 (dashed black lines),
location of the altimetry virtual station (circle, in black), and the
area used for extraction of flood extent (Lower Purus, pink polygons);
(b) Hydrological Response Units (Fan et al., 2015); (c) Bare Earth
Digital Elevation Model (O’Loughlin et al., 2016); (d) Köppen-Geiger
Climatic Zones (Kottek et al., 2006).
Hydrological-hydrodynamic model:
MGB
The MGB (“Modelo de Grandes Bacias”, a Portuguese acronym for “Large
Basin Model”) is a semi-distributed, hydrological-hydrodynamic model
(Collischonn et al., 2007; Pontes et al., 2017). It was chosen for this
study because (1) it has been
widely and successfully applied in several South American basins
(e.g., Paiva et al., 2013;
Siqueira et al., 2018); (2) it is representative and similar to other
conceptual hydrological models like VIC (Liang et al., 1994) and SWAT;
and (3) the hydrological component is tightly coupled to a hydrodynamic
routing scheme, allowing the simulation of complex flat, tropical
basins. Moreover, the source code of MGB is freely available at
www.ufrgs.br/lsh.
Within the model structure, basins are
discretized into unit-catchments,
which are further divided into Hydrological Response Units (HRU’s) based
on soil type and land use. A vertical water balance is performed for
each HRU, considering canopy interception, soil infiltration,
evapotranspiration, and generation of surface, subsurface and
groundwater flows. Soil is
represented as a bucket model with a single layer. Flow generated in
each HRU is routed to the outlet of the unit-catchment with linear
reservoirs. Outflow from each unit-catchment is then propagated through
the stream network by using a 1D hydrodynamic model based on the
inertial approximation proposed by Bates et al. (2010). The stream
network is derived from Digital Elevation Model (DEM) processing. The
model has 19 parameters, which are further detailed in the next section.
Other model inputs are precipitation, climate data, soil type and land
use maps, which are further described in section 2.6 Model Setup.
A priori uncertainty of model
parameters
Within MGB model, there are parameters related to vegetation cover
(albedo, leaf area index, vegetation height and Penman-Monteith surface
resistance), river hydraulics (Manning’s roughness, and width and depth
parameters related to geomorphological relationships), and conceptual
parameters related to soil water budget (Wm, b, Kbas, Kint, XL, CAP, Wc,
CI, CS, CB), which are further detailed in Supporting Information (Table
S2). Out of the 19 model parameters, six are fixed and 13 are
calibrated.
The a priori uncertainty of MGB model parameters is estimated based on
their variability as reported in literature (references in Table S2).
Supporting Information (Table S2)
presents the calibration parameters, their initial values, range, and
the references that support these assumptions.
Sensitivity analysis
In order to understand how different parameter sets (river hydraulic,
soil, vegetation) affect model output variables (river discharge, flood
extent, river water level, soil moisture, evapotranspiration and
terrestrial water storage), multiple model runs were conducted
considering four uncalibrated model setups: (1) varying only soil
parameters; (2) varying only vegetation parameters; (3) varying only
hydraulic parameters; (4) varying all parameters together.
One hundred runs were conducted in
triplicate to ensure that convergence is not dependent on the initial
parameter sets, thus resulting in 300 runs for each setup. In this step,
no RS-based calibration is performed yet.
Parameters were varied considering a uniform distribution, and results
were analyzed in terms of mean RMSD (root mean square deviation) of each
variable, by comparing each run with a reference one (i.e., the initial
run with the initial parameter set as defined in Supporting Information
Table S2). This was performed in order to understand the sources of
model uncertainties related to different sets of parameters (e.g., are
flood extent estimates sensitive to vegetation parameters, or are ET
estimates sensitive to hydraulic parameters?). The dispersion of model
outputs was also compared to uncertainty in the observations, as derived
from literature.
To understand which variables are inter-related in the model, we further
analyzed the results of setup “(4) varying all parameters together”.
This was done by firstly computing the Kling-Gupta Efficiency metric
(KGE; Gupta et al., (2009)) between the perturbed runs and a reference
one (i.e., run with the initial parameter set) for each variable, and
then calculating the Pearson correlation (r) between KGE values for each
pair of variables (e.g., discharge and water level, discharge and flood
extent, and so forth). This
experiment is relevant to evaluate whether two variables get improved or
get worsened together, or whether a variable improvement impacts on the
deterioration of another. In other words, this approach allows to
evaluate the correlation between the variables.
Model setup
The Bare Earth Digital Elevation Model (DEM; O’Loughlin et al., 2016)
(Figure 2c) was used for stream network computation and basin
discretization with the IPH-HydroTools GIS package (Siqueira et al.,
2016). The original DEM resolution
is 90 m, and it was resampled to 500 m to facilitate GIS processing. An
upstream area threshold of 100 km2 was adopted to
delineate the drainage network, and unit-catchments were discretized by
dividing the stream network into fixed reach length of 10 km, resulting
in 2957 unit-catchments for the whole basin. Soil type and land cover
maps were extracted from the HRU discretization developed by Fan et al.
(2015) (Figure 2b): (1) deep and (2) shallow forested areas, (3) deep
and (4) shallow agricultural areas, (5) deep and (6) shallow pasture,
(7) wetlands, (8) semi-impervious areas, and (9) open water, where
“deep soils” refer to soils with high water storage capacity, and
“shallow soils” are those with low water storage capacity. In the
Purus River Basin, 57.4% of the region is covered by forest with deep
soils, 26.9% by forest with shallow soils, and 13.7% by wetlands
(i.e., river floodplains). Daily precipitation data were derived from
TMPA 3B42 (version 7), with spatial resolution of 0.25º x 0.25º (Huffman
et al., 2007; available at:
<https://gpm.nasa.gov/data-access/downloads/trmm>),
and were interpolated with the nearest neighbor method for the centroid
of each unit-catchment. Long term
climate averages for mean surface air temperature, relative humidity,
insolation, wind speed and atmospheric pressure were obtained from the
Climatic Research Unit database (New et al., 2000; available at:
<http://www.cru.uea.ac.uk/data>), at a
spatial resolution of 10’, and also interpolated with the nearest
neighbor method.
Model calibration
The MOCOM-UA calibration algorithm (Yapo et al., 1998; Multi-objective
global optimization for hydrologic models) was adopted due to its
satisfactory performance when coupled with hydrological models (e.g.,
Collischonn et al., 2008; Maurer et al., 2009; Naz et al., 2014).
MOCOM-UA is an evolutionary algorithm, based on SCE-UA (Duan et al.,
1992), that simultaneously optimizes a model population with respect to
different objective functions. The
model population consists of randomly distributed points within the
parameter search space, and it reflects the a priori uncertainty of
model parameters. Here, the population size was set to 100 individuals.
The altered model parameters and their respective ranges are described
in Supporting Information Table S2. All calibration experiments are
repeated three times with differing initial parameter sets to ensure
that convergence is not dependent on the initial parameter sets.
Objective functions to be optimized depend on the calibration setup. In
the single-variable calibration, for each variable, three objective
functions (OF) that summarize the agreement between simulated
and observed (RS) time-series are simultaneously optimized: Pearson
correlation (r), ratio of averages
(\(\text{\ μ}_{\text{sim}}\ /\ \mu_{\text{obs}}\ \)), and ratio of
standard deviations
(\(\ \sigma_{\text{sim}}\ /\ \sigma_{\text{obs}}\ \)), which are
associated to the individual terms of KGE metric. These 3 objective
functions are depicted in Equations 1 to 3, where X denotes the assessed
variables (Q, h, A, TWS, ET or W).
\begin{equation}
\text{OF}_{1}=\ \left(\frac{\text{\ μ}_{\text{sim}}}{\mu_{\text{obs}}}\right)_{\text{X\ \ }}\left(1\right);\ \ \ \ \ \ \text{OF}_{2}=\ \left(\frac{\sigma_{\text{sim}}}{\sigma_{\text{obs}}}\right)_{\text{X\ \ }}\left(2\right);\ \ \ \text{\ \ OF}_{3}={\ r}_{\text{X\ \ }}(3)\nonumber \\
\end{equation}For the multi-variable calibration, the objective functions are the KGE
of each variable considered: firstly, five objective functions were
considered (KGE of all variables except discharge), as depicted in
Equations 4 to 8.
\begin{equation}
\text{OF}_{1}=\text{KGE}_{h}\ \left(4\right);\ \ \ \ \ \text{OF}_{2}=\text{KGE}_{A}\ \left(5\right);\text{\ \ \ \ OF}_{3}=\text{KGE}_{\text{TWS}}\ \left(6\right);\text{\ \ \ \ \ OF}_{4}=\text{KGE}_{\text{ET}}\ \left(7\right);\text{\ \ \ \ \ \ \ OF}_{5}=\text{KGE}_{W}(8)\nonumber \\
\end{equation}Secondly, two objective functions were adopted (KGE of selected
variables 1 (x) and 2 (y)), as depicted in Equations 9 and 10.
\begin{equation}
\text{OF}_{1}=\text{KGE}_{x}\ \left(9\right)\ \ \ \ \ \ \ \ \ ;\text{\ \ \ \ \ \ \ \ \ \ OF}_{2}=\ \text{KGE}_{y}(10)\nonumber \\
\end{equation}Results are expressed in terms of a Skill Score (S) (Equation 11; Zajac
et al., 2017), in order to evaluate the improvement (or
deterioration) in the representation of a variable when the model is
calibrated with a given variable, compared to the uncalibrated setup.
\begin{equation}
S=\frac{\text{KGE}_{\text{calibrated}}-\text{KGE}_{\text{initial}}}{1-\text{KGE}_{\text{initial}}}\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ (11)\nonumber \\
\end{equation}KGEcalibrated is the mean KGE resulting from running the
model with the calibrated parameters. KGEinitial is the
mean KGE resulting from running the model with the a priori parameter
sets.
Calibration/Evaluation Data
In the next paragraphs we introduce the data used for model calibration
and evaluation, as well as how MGB outputs were evaluated in comparison
to them.
-In-situ discharge measurements were obtained from the Brazilian
Water Agency Hidroweb database (available at:
<http://www.snirh.gov.br/hidroweb/publico/apresentacao.jsf>),
at the gauge “Canutama” (code: 13880000; location: S ° 32’ 20.04”; W
64° 23’ 8.88”; drainage area: 236,000 km², period of data availability:
1973 to 2016). Uncertainty in discharge observations can be estimated as
ranging from 6.2% to 42.8% at the 95% confidence level, with an
average of 25.6% (Di Baldassarre & Montanari, 2009). Discharge was
evaluated on a daily basis.
- Remotely sensed water level data were obtained from Jason-2
mission, which presents an orbit cycle of approximately 10 days, and
tracks separated by approximately 300 km at the equator (Lambin et al.,
2010). It presents an accuracy of approximately 0.28 m (Jarihani et al.,
2013), and data are available since 2008. The virtual station presented
in Figure 1 corresponds to Track number 165. Processed data for this
study were downloaded from the Hydroweb/Theia database (available at:
<http://hydroweb.theia-land.fr>). Water level was
computed in MGB at the unit-catchment associated to the altimetry
virtual station, being an advantage of using the hydrodynamic scheme for
flood routing instead of the Muskingum simplification. Simulated and RS
water level data were compared every 10 days in terms of anomaly (values
subtracted from long term average).
- Satellite flood extent
data were derived from ALOS-PALSAR imagery, which presents a ground
resolution of 100 m (Rosenqvist et al., 2007). Images were downloaded
from Alaska Satellite Facility (available at:
<https://www.asf.alaska.edu/>) in
processing level 1.5, which already presents geometric and radiometric
corrections. A 3 x 3 median filter was used to remove speckle noise (Lee
et al., 2014). Images were classified into water (backscattering
coefficient less than -14 dB), non-flooded forest (between -14 dB and
-6.5 dB), and flooded forest (higher than -6.5 dB) classes, according to
Hess et al. (2003) and Lee et al. (2014).
The uncertainty of flood extent
estimates was estimated based on the RMSE between the resulting
classification of this study, and the dual-season mapping developed by
Hess et al. (2003). Simulated and RS flood extent data were compared for
the pink area depicted in Figure 1, in order to avoid spurious flood
extent data in regions that are known to be not subject to flooding.
ALOS-PALSAR presents a recurrence
cycle of 46 days (from 2006 to 2011), and flood extent data were
available and compared to MGB for 21 dates.
- Satellite-based terrestrial water storage (TWS) anomalies were
extracted from GRACE mission, launched in March 2002. GRACE provides
monthly TWS estimates based on anomalies in gravitational potential, at
a resolution of 300-400 km, with a uniform accuracy of 2 cm over the
land and ocean regions (Tapley et
al., 2004). TWS anomalies were retrieved from three processing centers -
GFZ (Geoforschungs Zentrum Potsdam, Germany), CSR (Center for Space
Research at University of Texas, USE), and JPL (Jet Propulsion
Laboratory, USA), available at
<https://grace.jpl.nasa.gov/>,
and then the mean value based on
the three products was averaged for the whole basin. In MGB, TWS values
were computed as the sum of water storage of all hydrological
compartments: river, floodplains, soil, groundwater and vegetation
canopy. Simulated and RS-based TWS were compared in terms of anomaly
(values subtracted from long term average) at a monthly time-scale.
- Satellite-based evapotranspiration estimates were retrieved
from the MOD16 product, derived by an algorithm presented by Mu et al.
(2011) based on the Penman-Monteith equation. The dataset covers the
period 2000-2010 with a spatial resolution of 1 km for global vegetated
land areas. Because of that, even though MGB evapotranspiration is
calculated for flooded areas (open water evaporation in main channel and
floodplains) and vegetation for water balance purposes, only the
vegetation-ET output was compared to MOD16. MOD16 products are provided
in 8-days, monthly and annual intervals. Monthly intervals were used
here and averaged for the whole basin. Accuracy of MOD16 along the
Amazon basin is estimated as 0.76 mm/day (Gomis-Cebolla et al., 2019).
MOD16 data is available at: <
https://www.ntsg.umt.edu/project/modis/mod16.php>. In
MGB, evapotranspiration is computed via Penman-Monteith equation, based
on the climate input variables.
- Satellite-based soil moisture is derived from the SMOS mission
(Kerr et al., 2001), processed by the Centre Aval de Traitement des
Données SMOS (CATDS), and downloaded in processing level 4, which
combines lower level products with data from other sensors and
modeling/data assimilation techniques. The daily L4 root zone soil
moisture product at 0-1 m soil depth (Al Bitar et al., 2013) were used
(available at:
<https://www.catds.fr/Products/Available-products-from-CEC-SM/L4-Land-research-products>),
and data from ascending and descending orbits were averaged for the
whole basin. In MGB, soil moisture as a saturation degree was computed
as the water in the soil compartment divided by the maximum water
capacity of the soil (Wm parameter). Since MGB estimates
saturation degree values for a soil bucket reservoir, SMOS values were
rescaled for the range 0 - 100% according to the Min/Max Correction
method described by Tarpanelli et al. (2013) and applied by some studies
(e.g., Rajib et al., 2016; Silvestro et al., 2015), and them compared to
MGB on a daily time-scale as an average for the whole basin.