INTRODUCTION

Wetlands are saturated with water, either permanently or seasonally and present structural and functional properties different from other ecosystems. They regulate water and biogeochemical cycles and provide several ecological services like water purification, flood control and biodiversity preservation, among others (Bucher, 1999; Bedford et al., 2001; Keddy, 2010). Wetlands are the more altered ecosystem on Earth due to human activities (Costanza et al., 1997; Mitsch and Gosselink, 2000; Dudgeon et al. 2006). In many cases, the productive use of wetlands is based on its transformation to terrestrial ecosystems (e.g. embankments), which alters drastically their hydrological functio­ning and thus its capacity to provide ecological services (Bunn and Arthington, 2002; Mitsch, 2014).

The Delta of the Paraná River surrounds the most popu­lated area of Argentina. More than fifteen million people benefit from its ecological services such as water supply and recreation, then, environmental conflicts product of its management are clearly visible by society. Extensive far­ming is the most traditional productive activity. However, in the last two decades, livestock has intensified and new activities appeared such as forestry, agriculture or urba­nization, among others. These new activities are often accompanied by infrastructure projects that alter the hydrodynamics of the system (Donadille et al., 2010). Even though wetlands are resilient systems (Quintana et al., 2008), the new activi­ties already caused environmental problems (as the great fires of 2008; Sione et al., 2008). It is necessary to define levels of human intervention that preserve the ecological integrity of the systems for which it is necessary to have information about its hydrodynamic process.

Similar to other wetlands of the world, structure and functioning of the Paraná Delta varies widely in space and time (Bunn and Arthington, 2002). Alternating periods of high and low water and occurrence of extraordinary flooding events are part of the natural dynamics that help to maintain their ecolo­gical integrity (Bó and Malvárez, 1999). In turn, the interaction between hydrological and geomorphological processes defines a wide variety of environments with its own functional characteristics (Malvárez, 1999). Currently, there are regional flooding forecasting systems (provided by the Water National Institute) and static maps of landscape and vegetation units (Malvárez, 1999; Salvia et al., 2007, Morandeira et al., 2011). However, to deal with such diverse ecosystems, a tool that provides timely information on which converge detailed spatial and temporal resolutions is needed. In that sense, there are studies of vegetation and inundation dynamics based in NDVI time series (Salvia, 2010; Borro et al., 2014).

Remote sensors allow carrying out studies of large a­reas at detailed spatial and temporal scales. There are some products specifically developed to register wetlands flooding. Some of these products are the Small Water Bodies monitoring product (Bartholomé, 2008) and the NASA Earth System Data Record of inundated wetlands (NASA JPL wetland, 2013), among others. However, none of these systems meet simultaneously the three desired attributes required for this kind of ecosystems: fine temporal and spatial scales and near real time availability. Papers describing the regional hydrodynamics agree that the Moderate Resolution Imaging Spectroradiometer (MODIS) sensor, on board on AQUA/TERRA platforms, is the most suitable for this kind of studies (Brakenridge and Anderson, 2006; Sakamoto et al., 2007; Islam et al., 2009; Chen et al., 2014; Pekel et al., 2014). Some models based on different spectral indices derived from MODIS have been proposed (Sakamoto et al., 2007; Chen et al., 2014; Pekel et al.; 2014). These indices rely on the absorbance of water at the near infrared (NIR) and shortwave infrared (SWIR) bands of the electromagnetic spectrum (Brakenridge and Anderson, 2006). As some components of water that vary spatially and temporally (e.g. phytoplankton, or dissolved organic matter) affect spectral properties (Arst, 2003; Pekel et al., 2014) it is necessary to assess which indices are best suited for the study region, where water and sediment pulses can influence the turbidity and spectral properties of water (Minotti and Borro, 2011).

The objective of this work was to develop a water cover estimation model based on MODIS information to describe the inundation patterns of the Paraná Delta. First, we evaluated the accuracy of different vegetation indices to estimate water cover under different water scenarios. Then, the best model was used to describe the inundation patterns for a period of 12 years every 16 days. This provides a tool towards a more sustainable production in the Paraná Delta Region and a base reference for other similar areas of the world.

MATERIALS AND METHODS

2.1 Study area

The Paraná Delta covers an area of around 17500 km² (Malvárez, 1997; figure 1). It is the terminal portion of the Paraná River Basin, which is the second largest in South America (2.6 Mkm²). The hydrology of the Paraná Delta is mainly influenced by the hydrological regime of the Paraná River, which depends on precipitations that are higher during summer months. The hydrology of the Paraná Delta is also influenced by the hydrological regimes of the Uruguay and Gualeguay rivers and by astronomical or meteorological tides of the Río de la Plata River, which affect di­ffe­rent sectors (Mujica, 1979, Quintana and Bó, 2011). Local rainfall is less important, and affects only those sectors not influenced by flooding from rivers (Quintana and Bó, 2011).

Figure 1Map of the Paraná Delta River Region, with the landscape units defined by Malvárez (1999, indicated with letters and numbers). Blue arrows indicate the main rivers that affect the hydro-patterns of the Paraná Delta (their size is proportional to the flow of the rivers). The city of Victoria is also marked. Top-right map: location of the study region (rectangle) within Argentina (dark gray) and the rest of South America (light gray). The Paraná River is marked with a blue line and its watershed is shadowed. 

At regional scale, different landscape units are associa­ted to pography and hydrological regime (Malvárez, 1997; figure 1). Landscape units share a common heterogeneity pattern with uplands at the perimeter of the islands, medium lands and frequently flooded lowlands, although their proportion varies between units. The unit A is dominated by non-vegetated lagoons surrounded by semi-permanently flooded and meander plains. The unit B is dominated by NO-SE rivers surrounded by uplands, and lagoons su­rrounded by marshes. The unit C is dominated by parallel sand ridges separated by depressions into three units: ridges with grasslands (C1), ridges with forests (C2), and ridges with forests and streams that drain to the Paraná River (C3). The unit D is dominated by a very flat semi-permanently flooded plain with lagoons and streams. The unit E presents uplands due to meanders, whose turns enclose permanent lagoons. The unit F is an extensive and low drained flood-plain formed during the Holocene marine ingression. The unit G is formed by ancient deltas of the Nogoyá, Clé and Gualeguay rivers. The H unit is formed by the Ibicuy Island that was included in the region after sea regression. The lower delta unit (which includes the Paraná and Uruguay River units) is the only strictly deltaic region, with a bidirectional hydrological regime. A great proportion of this unit was drained for forestation.

2.2 Water cover estimation

As ground truth we used 30 m spatial resolution water presence-absence maps developed by Sepulcri et al. (2012) from eight LANDSAT images from different water scenarios defined from the annual average flow of the Paraná River: low water: 10/01/2009, 02/05/2009, normal: 23/02/2003, 25/05/2003, 04/03/2005, 05/04/2005 and high water: 29/05/2007, 12/01/2010. These maps were calculated by photointerpretation following a non-supervised “isodata” classification (maximum 30 classes, Ball and Hall 1965) based on three normalized vegetation indices derived from LANDSAT images: NDVI, NDSI and NDWI, based on the 4 (750-900 nm) and 3 (630-690 nm), 5 (1550-1750 nm) and 4, and 2 (525-605 nm) and 7 (2090-2350 nm) LANDSAT bands, respectively. To assess accuracy of LANDSAT based maps we compared them with ten 3 x 3 km² high resolution images from close dates (less than 10 days) and where water height remained constant. Images were extracted from Google Earth Pro and classified into water presence-absence by photointerpretation. Conside­ring Google Earth images as ground truth, mean (± standard deviation) underestimation and overestimation errors of LANDSAT based maps were 5.2% (±4.5) and 2% (±1.5) respectively. For image processing we used the GRASS software (http://grass.osgeo.org/).

To merge the ground truth (LANDSAT water presence-absence maps) and MODIS images we overlapped the 250 m MODIS grid with the water presence-absence maps and calculated the proportion of LANDSAT pixels classified as water inside each MODIS pixel. To avoid confounding temporal and spatial variation we selected only those pixels that were within 6 water cover ranges (0-10, 10-30, 30-50, 50-70, 70-90, 90-100%) for the 8 dates. They covered 24% of the area (21% corresponded to the range 0-10, 2.5% to the range 90-100 and 0.5% corresponded to the 4 intermediate ranges). We evaluated as response variables either water cover inside MODIS pixels (sub-pixel resolution) or water presence-absence (pixel resolution). To define water presence-absence we evaluated different pixel water proportion thresholds.

As water cover or presence estimators we evaluated several water related spectral indices previously used to assess water, and the color components expressed in the HSV space (Hue, Saturation and Value) derived from MODIS bands (table 1). According to the HSV model, H and S represent the chromaticity, while V represents the brightness. H is the perceived color and is expressed in degrees, where 0º or 360º represent the red, 120º the green, and 240º the blue colors. In turn, S express how far from the grey is a color at a given V (Pekel et al., 2014). These indices are mainly based on the liquid water absorption features in the infrared range (around 970, 1200, 1450 and 1950 nm, Cu­rran 1989, Asner 1998). We also evaluated an index based on temperature (i.e. the diurnal land surface temperature difference, based on the MODIS 11A2 product, Wan et al., 2002), but due to it lower spatial resolution (1 km), and not better performance, we decided not to show the results.

Table 1Spectral indices and HSV models utilized for the water cover estimation models. NIR is the surface reflectance value in the near infrared (841-876 nm, MODIS Band2), Red (620-670 nm, Band 1), Blue (459-479 nm, Band 3), NIR2 (1230-1250 nm, Band 5), and SWIR is the short-wave infrared (1628-1652 nm, Band 6), and MIR is the medium infrared (2105–2155 nm, Band 7). 

To calculate water cover or presence estimators, we obtained the MODIS pixel values for the required bands and for all the studied area from the MODIS 16 days composite that included the eight LANDSAT images defined above. The 250 m spatial resolution bands (Red and NIR) were obtained from the MOD13Q1 16 days composite, while the 500 m spatial resolution bands (Blue, NIR, NIR2 and SWIR) were obtained from the MOD09A1 8 days compo­site and were resampled to 250 m with the nearest neighbor method. All the MODIS pixels for the eight dates presented “Good Data” according to the “pixel reliability” layer, which contains ranked values describing overall pixel quality.

To evaluate the best water cover or presence estimator we developed regression and tree models with the R software (R Core Team 2013). For water cover estimation (sub-pixel resolution), regression models included simple or multiple estimators and linear or non-linear relationships. As there can be a saturation response between water cover and the spectral indices used as estimators (Tucker, 1977) linear models included thresholds of the estimators. The optimal number of predictors for multiple regression models was assessed by the Bayesian Information Criterion (BIC. Schartz, 1978) using the “regsubsets” function of the leaps “R” package (Lumley 2017). The BIC penalizes models by the number of parameters to reduce over-fitting. Non-linear logistic models were performed with the “glm” function of the MASS “R” package (Venables and Ripley, 2002). Tree models were performed with the “tree” function of the Party “R” package (Hothorn et al., 2006). They included all the estimators with 2, 4 and 8 final nodes. We did not include a higher number of final nodes because they lead to over-fitting errors. For water presence-absence estimation (pixel resolution), we used tree models. To evaluate accuracy of both water cover and water presence-absence models we used a re-sampling method. First, we built the above mentioned models based on 50 pixels randomly selected from each of four water cover ranges (0-25, 25-50, 50-75 and 75-100%) and date (N = 1600 pixels). Second, we tested those models with 1600 pixels randomly selected into the whole image. We used the coefficient of determination R² to evaluate model fit, and root mean square error (RMSE) and Kappa index to evaluate model accuracy. We recorded the values of the evaluators for each resampling (N = 1000) and finally, we calculated mean and standard deviation va­lues to compare and select the best model.

2.3 Inundation patterns assessment

After choosing the best model, we characterized the inundation patterns of the study area. The study period extended from January 2001 to December 2012 and integrated 276 16-days MODIS images. We only included data for pixels with “Good Data” according to the pixel reliability band. Besides, 16 ima­ges with a high proportion (>30%) of low quality pixels were completely excluded; 3 images for 2012 and between 0 and 2 for the rest of the years. During the study period there were two extraordinary flooding events that exceeded evacuation levels on the riparian cities. The 2007 event extended from March to September (12 MODIS scenes), while the 2010 event exten­ded from July 2009 to September 2010 (25 MODIS scenes).

With the remaining 260 “good quality” MODIS scenes, we generated water presence-absence maps to characterize the inundation patterns by two ways. The first way consis­ted on overlapping water presence-absence maps with the landscape units’ map (Malvárez, 1999; figure 1). Then, we described average water cover and intra- and inter-annual variation for each landscape unit. The second way consis­ted on developing a classification based on the water pre­sence-absence maps and relating water cover with river height. The classification was done in three steps. First step, we obtained the two extreme classes which included those pixels never classified as water and those classified as water most of the dates (> 80% of the MODIS scenes). Second step, we obtained three classes that included pixels classified as water only during the extraordinary flooding events: only during 2007, only during 2010 or in both events. Third step, we developed an unsupervised classification based on water pre­sence for the remaining area. We performed a principal component analysis based on water presence of the 223 MODIS scenes that did not include extraordinary flooding events to reduce data dimensionality (GRASS i.pca command). Then, with the firsts three axis, that explained 43% of the variation, we obtained five classes by an “isodata” clustering algorithm (GRASS i.cluster and i.maxlik commands). Once obtained the classification, we evaluated the relationship between water height of the Paraná river at the location of Victoria (figure 1) and the percentage of pixels classified as water for a 7x7 moving window as a first approach to relate water cover with a variable related to hydrodynamics process.

RESULTS AND DISCUSSION

3.1 Water cover estimation based in MODIS

Water cover into MODIS pixels was associated to spectral indices (figure 2). As described by Guershman et al. (2011), NDVI and EVI presented a non-linear negative relationship with water cover inside pixels (figures 2 A-B). It was observed a high variability between dates, especially at low water cover values, which could be associated to the sensibility of both indices to changes in the vegetation cover and status (Baret and Guyot, 1991; Gamon et al., 1995; Van Leeuwen and Huete, 1996; Huete et al., 2002), which have a higher weight when water cover is low. This pattern shows the limitations of both indices for the objective of this work (to estimate water cover), but also warns about their limitation to estimate aspects of vegetation (e.g. fraction of incoming radiation absorbed by the vegetation, primary productivity) on systems of highly variable water cover. In these systems, water level and the amount of water above and below water highly affect the spectral signal (Beget et al., 2013) and other spectral indices can be more accurate (i.e. the Visible Atmos­pheric Resistant Index, Gosh et al., 2016).

Figure 2Relationship between water cover inside MODIS pixels derived from LANDSAT images and MODIS spectral indices. Each point represents the average value of all the pixels inside each water cover range (0-10, 11-30, 31-50, 51-70, 71-90, 91-100%) and date. Legend in the top left indicates MODIS composite dates containing LANDSAT images; each date is represented by a different symbol. 

From the three indices previously used to detect water (figure 2 C-E), the NDWI1 was the most closely related with water and showed an asymptotic relationship. The NDW1 was less variable between dates than the NDVI and the EVI, probably associated to the higher sensitivity to water cover (Beget et al., 2013) and lower sensitivity to vegetation status of the shortwave infrared compared to the near infrared (Curran, 1989; Asner, 1998).

Between the three components that define the color (fi­gure 2 F-H), the H index (associated to the perceived color) was the closely related to water cover. As NDWI1, the relationship between water cover and H was asymptotic. However, the range in which H remained sensitive to water cover was more limited, which reduced the capacity of H to discriminate intermediate water cover levels. At maximum water cover, H was around 160º which means higher re­lative values of near infrared and red bands compared to medium infrared and coincides with the higher water absorption in the medium infrared (Beget et al., 2013). The V associated to brightness captured better the difference between dates than between water cover ranges. This could be due to factors that affect the complete scene like sun and observation angles, or atmospheric effects rather that land cover type (e.g. Pekel et al., 2014).

The capacity to predict water cover of the linear with thresholds, non-linear, and tree models was similar to that found in previous works (Huang et al., 2014) although higher errors were observed at intermediate water cover levels (figure 3, table 2). Predictive accuracy was not improved by increasing the complexity of models. Including other indices in the regression models based on NDWI1, or more than four final nodes in tree models did not improve model accuracy significantly (p > 0.5). The coefficient of determination (R²) and the prediction error (RMSE) between observed and predicted values was not different (p > 0.01) between models (table 2). However, the Kappa coefficient was slightly higher (p=0.09) for the sigmoid than for the tree model, but was not different from the linear with thresholds model. This difference was due to the higher number of pixels in the extreme water cover ranges where sigmoid and linear with thresholds models were more accurate.

Figure 3A) Relationship between water cover inside MODIS pixels and NDWI1 for 20000 randomly selected pixels. Darker colors indicate a higher density of pixels. Dashed lines represent the best fit lines for linear (Water Cover = 66 x NDWI1 + 57) and sigmoid (Water Cover = [exp(0.86 + 4.6 x NDWI1)] / [1 + exp(0.86 + 4.6 x NDWI1)]) models. B) Relationship between H (left) and NDVI (right) with NDWI1 for 20000 randomly selected MODIS pixels. Color intensity indicates water cover inside pixels. Partition of the x axis and horizontal lines indicate the threshold values of the tree model separating water cover levels (first node: NDWI1 = -0.29, second node: H = 119 and third node: NDVI = 0.37). 

Water presence or absence estimation models, even though resigned sub-pixel resolution, improved the predictive accuracy compared to water cover estimation models. Due to the higher simplicity and accuracy of the water presence-absence estimation models we considered that they were more adequate to assess regional inundation patterns (see section 3.2). The red and SWIR bands, combined in the NDWI1, capture most of the information of water. The predictive accuracy of models based in a NDWI1 threshold was not significantly improved by models including a higher number of estimators and was similar to that obtained in other studies using more complex decision rules (Sakamoto et al., 2007; Ordoyne et al., 2008; Guershman et al., 2011; Pekel et al., 2014).

Table 2Water cover classification accuracy for the three kind of models: linear with lower an upper threshold, sigmoid and tree with four final nodes. Accuracy was assessed by the coefficient of determination (R2), the prediction error expressed by the root mean square error (RMSE) and the Kappa coefficient (k). To calculate k, water inside MODIS pixels was aggregated into four ranges (0-0.25, 0.25-0.50, 0.50-0.75 and 0.75-1.0). Values represent the mean (and standard deviation) for 1000 re-samples. 

Water cover threshold to classify a MODIS pixel as water or non-water affected predictive accuracy and bias. Underestimation and overestimation errors decreased as the threshold increased, but the overestimation was higher at low threshold values and decreased more steeply than the underestimation error (figure 4). The Kappa index increased monotonically until threshold values around 60% (which coincides with a NDWI1 of -0.2), where underestimation and overestimation errors were equilibrated. An optimal threshold of water cover between 50% and 80% was also found in other wetland region (Ticehurst et al., 2015). Besides being equilibrated in underestimation and overestimation errors, a 60% of the water cover is enough to affect the possible land management decisions. In forest with 40% of tree co­ver, a 60% of water cover would mean that the understory is completely covered with water.

Figure 4Changes in prediction error (black lines), NDWI1 threshold (gray line) and Kappa coefficient (dotted gray line) for water pre­sence or absence estimation models based on NDWI1 according to the water cover inside a pixel MODIS used to define the pixel as water. 

3.2 Inundation patterns of the Paraná Delta River

The estimation of water presence-absence based on NDWI1 every 16 days for 12 years showed differences between landscape units (Malvárez, 1999) on mean water cover and on its temporal variation (figure 5). We aggregated landscape units into three groups according to inundation patterns: 1) those with high within year water co­ver variation (standard deviation > 2.5%, units A, B, F and G, figure 5 A), 2) those with high within year water cover variation mainly in years with extraordinary flooding (water cover during flooding more than 5 times higher than during normal years, units D, E, H and Paraná de las Palmas N., figure 5 B), and 3) those having low water cover variation within and among years Delta (Paraná de las Palmas S., Paraná Guazú and Uruguay River, figure 5 C). Units of the first group were located at the north-east portion. From these, the unit B where temporal lagoons are more frequent (Malvárez, 1999; Borro et al., 2011) was the only one that decreased water cover during the 2008 drought. Units of the second group were located at the south-east. The high variation during flooding years is particularly evident in the C3 unit, which is the only that drains to the Uruguay River. Units of the third group were located at the Lower Delta. They are influenced by the Paraná and Río de La Plata Rivers. Extraordinary flooding events were not evident on these units. Instead, a flooding event of the Uruguay River that occurred during the 2009 is evident in the Uruguay Ri­ver unit. The proportion of different vegetation units, which are distributed and highly fragmented across the region (Salvia et al., 2007) could explain this inundation patterns.

Figure 5Mean +/- intra-annual standard deviation for water cover for each year (2001-2012) and vegetation unit (Malvárez, 1999). Years with extreme flooding events (2007 and 2010) are marked in thick gray line. Vegetation units are divided into three plots according to its location: A) north-east (units A, B, F and G), B) south-east (units C, D, E, H and Paraná de las Palmas N.) and C) lower Delta (Paraná de las Palmas S, Paraná Guazú and Río Uruguay). 

The land cover classification based on water presence or absence maps showed the strong variability of the study area (figures 6 and 7). Extreme classes (figure 6 A) co­vered only 27% of the area; 21% was never covered by water while 6% was covered by water more than 80% of the dates. The 3 classes that only included water during the extraordinary flooding events covered 14% of the area (fi­gure 6 B). In contrast, most of the area (58%) was covered by pixels that presented water not only during extraordinary flooding events, which were aggregated on five classes (fi­gure 7). Classes that presented water a few dates moreover extraordinary flooding events (classes A and B) co­vered the larger area. The Class A covered mainly the Lower Delta and the right margin of the Paraná River, while the Class B covered the central zone of the Medium Delta. The Class B showed higher sensitivity to flooding events than the Class A, expressed by an earlier incidence of the 2010 flooding (that begins at the end of 2009), and a larger water cover during 2011. Classes C, D and E presented a larger mean water cover. Classes C and D were located at the north of the Lower and Medium Paraná, respectively, while the Class E was spread across all the region, associated to water courses. On class C, there were several peaks similar to the extraordinary flooding events, but there were also periods were water cover was near to zero. In contrast, on Class D, although both extraordinary flooding events were more notorious, water cover was less frequently near to zero. Class E presented the higher mean water cover. Always at least 30% of the pixels were covered by water and extraordinary flooding events were not notorious.

Figure 6Maps showing the two extreme classes (A) and the three classes representing areas covered with water only during extreme flooding events obtained (B). The area never covered with water and more than 80% of the dates was 4695, 1632 km², while the area covered with water only during 2007, 2010 and during both years was 1486, 480 and 1406 km², respectively. 

Figure 7Left) Maps showing the five classes that varied its water cover not only during extraordinary flooding events obtained in the classification based on water cover maps. Right) Water cover along the year for each of the twelve years (2001-2012) and class (A, B, C, D and E). Years with extreme flooding events are marked with tick lines (2007 full and 2010 dashed). The area covered by the classes A, B, C, D and E was 6352, 3996, 1067, 1061 and 696 km² respectively. 

Some studies showed an association between time series of water area and water level in lakes (Zhu et al., 2014; Hu et al., 2015) and marshlands (Ordoyne et al., 2008). Similarly, we found a correlation between water level of the Paraná River at Victoria and water cover (figure 8). Particularly for the classes D, B and A from this work (figure 7), water level at Victoria explained, in average, 52%, 47%, and 30% of the variations in water cover respectively. This analysis allowed to identify those areas more linked to the Paraná River level fluctuations.

Figure 8Relationship between time series of water height of the Paraná river at the location of Victoria and the percentage of pixels classified as water for a 7x7 pixels moving window. Red saturation indicates the coefficient of determination (R2) of linear models between water height and water cover of each 7x7 pixels window. Areas never covered by water and more than 80% of the dates are marked in gray and black respectively. 

CONCLUSIONS

A tool that provides frequent and spatially explicit water cover estimations was developed and used. Many water related spectral indices were evaluated. Even though none of them accurately estimated water cover inside pixel MODIS (sub-pixel resolution), water presence or absence (pixel resolution) was accurately estimated by a model based on a MODIS NDWI1 threshold, which accurately separated pixels with less and more than 60% of water cover. More complex models did not improve accuracy.

By using the model, it was possible to map water cover in the region every 16 days for 12 years at 5 ha spatial resolution. Mean water cover and its temporal variation of previously defined landscape units (Malvárez, 1999) was assessed. Besides, a new classification based on the frequency and incidence of flooding events was developed and related to water level. These results demonstrated the utility of MODIS time series to characterize the hydrological patterns of the Paraná Delta. Further, the availability of new high- spatial and -temporal resolution radar and optical sensors (i.e. the Sentinel mission) can help to describe dynamics at a finer grain scale.

Spatially explicit water cover estimation can complement existing static information and in situ and monitoring stations. For example, the relationship of water cover with water level at different locations could improve flooding alert systems. Also, the combination of water cover and forage productivity estimations could allow having an idea of the temporal variation in the regional carrying capacity due to changes in the usable area. Besides, as this information can be delivered at almost real time, it could help for decision making during extraordinary flooding events.