Applying Independent Component Analysis on Sentinel-2 Imagery to Characterize Geomorphological Responses to an Extreme Flood Event near the Non-Vegetated Río Colorado Terminus, Salar de Uyuni, Bolivia

In some internally-draining dryland basins, ephemeral river systems terminate at the margins of playas. Extreme floods can exert significant geomorphological impacts on the lower reaches of these river systems and the playas, including causing changes to flood extent, channel-floodplain morphology, and sediment dispersal. However, the characterization of these impacts using remote sensing approaches has been challenging owing to variable vegetation and cloud cover, as well as the commonly limited spatial and temporal resolution of data. Here, we use Sentinel-2 Multispectral Instrument (MSI) data to investigate the flood extent, flood patterns and channel-floodplain morphodynamics resulting from an extreme flood near the non-vegetated terminus of the Rio Colorado, located at the margins of the world’s largest playa (Salar de Uyuni, Bolivia). Daily maximum precipitation frequency analysis based on a 42-year record of daily precipitation data (1976 through 2017) indicates that an approximately 40-year precipitation event (40.7 mm) occurred on 6 January 2017, and this was associated with an extreme flood. Sentinel-2 data acquired after this extreme flood were used to separate water bodies and land, first by using modified normalized difference water index (MNDWI), and then by subsequently applying independent component analysis (ICA) on the land section of the combined pre- and post-flood images to extract flooding areas. The area around the Rio Colorado terminus system was classified into three categories: water bodies, wet land, and dry land. The results are in agreement with visual assessment, with an overall accuracy of 96% and Kappa of 0.9 for water-land classification and an overall accuracy of 83% and Kappa of 0.65 for dry land-wet land classification. The flood extent mapping revealed preferential overbank flow paths on the floodplain, which were closely related to geomorphological changes. Changes included the formation and enlargement of crevasse splays, channel avulsion, and the development of erosion cells (floodplain scour-transport-fill features). These changes were visualized by Sentinel-2 images along with WorldView satellite images. In particular, flooding enlarged existing crevasse splays and formed new ones, while channel avulsion occurred near the river’s terminus. Greater overbank flow on the floodplain led to rapid erosion cell development, with changes to channelized sections occurring as a result of adjustments in flow sources and intensity combined with the lack of vegetation on the fine-grained (predominantly silt, clay) sediments. This study has demonstrated how ICA can be implemented on Sentinel-2 imagery to characterize the impact of extreme floods on the lower Rio Colorado, and the method has potential application in similar contexts in many other drylands.


Introduction
In some internally-draining (endorheic) dryland basins, ephemeral river systems decrease in size downstream and terminate at the margins of playas [1][2][3]. During rare flood events, flow and sediment may spill beyond the confines of the channel banks, and spread across extensive, fine-grained floodplains and the commonly saline playa surfaces [4,5]. These lower river reaches and the playas have received increasing research attention owing to their ecological, hydrological and geomorphological significance [6][7][8][9][10]. Some terminal dryland river and playa sedimentary systems may also provide modern analogs for unconventional hydrocarbon reservoirs [11,12]. Landsat time-series and high resolution satellite imagery (WorldView, QuickBird) have revealed the dynamic nature of fluvial landforms in these settings, which may include active and abandoned channels, crevasse channels and splays, and erosion cells (floodplain scour-transport-fill (S-T-F) features) [13,14]. This dynamism has short-term effects on flooding and sediment dispersal patterns, as well as longer term effects on the fluvial landscape and sedimentary successions [14][15][16].
In remote playas, rare flood events may be beneficial for ecohydrology, but in drylands more generally as well as in other hydroclimatic settings, floods also can be one of the most significant and widespread hazards [17]. In particular, large or extreme floods disrupt riparian ecosystems and present great threats to human life, properties, and other infrastructure, and hydrological and hydraulic models are needed to predict their impacts. Modeling of flood propagation and extent, both in natural and anthropogenically-induced events such as dam-break/dam-breach floods, requires accurate calibration and validation datasets [18][19][20][21]. However, there may be practical difficulties in acquiring the requisite field data on flood characteristics (e.g., flow duration, extent, distribution) to enable full model assessment. This may be particularly the case in low-gradient drylands where floods can be of great spatial extent and significantly restrict field access.
Consequently, remote sensing techniques, particularly those using optical imagery, have proven useful for mapping flood extent and flood patterns efficiently in terms of time and cost. Moisture sensitivity on visible, near infrared (NIR) and short-wave infrared (SWIR) bands enables land and water to be distinguished. In addition, the daily data acquisition of coarse resolution sensors (i.e., greater than 250 m) such as MODIS, VIIRS, and OLCI/SLSTR makes it possible to map flooding immediately following a rainfall-runoff event [22,23], especially as data may be collected at different times of the day. However, variable vegetation and cloud cover, as well as spatial and temporal data resolution issues may undermine the accuracy of results. Furthermore, change detection methods such as supervised classifiers require a priori knowledge for target classes to create training data, which is time consuming and potentially introduces human errors in pre-processing. Recently Chignell et al. demonstrated the usefulness of the combination of Landsat Operational Land Imager (OLI) data and independent component analysis (ICA), an unsupervised image classifier, for flood mapping [24]. Originally derived from the fields of signal processing and neural computing, the statistical technique ICA decomposes data into statistically independent sources with no a priori knowledge [25,26]. Chignell et al. mentioned the strengths of ICA such as concomitant delineation of flood-related surface water and soil, and also separation of non-flood change classes [24]. However, their study found Remote Sens. 2018, 10, 725 3 of 19 difficulties in delineating flood edges due to the low spatial resolution of the OLI imagery (30 m). High cloud cover immediately after floods may also undermine the accuracy of flood extent mapping undertaken using Landsat and other optical imagery; selection of post-flood, cloud-free satellite imagery may be possible, although for some applications it might be important to use imagery as close to the flood peak as possible.
Recently available Sentinel-2 Multispectral Instrument (MSI) imagery features relatively high resolution in spectral bands of visible to NIR (10 m) and red edge to SWIR (20 m) [27]. Owing to the improved spatial resolution, visible, NIR and SWIR bands have been proven useful for accurately extracting flooding extent (e.g., [24,28]). Modified normalized difference water index (MNDWI) using Sentinel-2 images has been successfully applied to extract surface water bodies [28,29]. By contrast, less attention has been directed to the impacts of single peak flood events on channel-floodplain geomorphology, despite the potential influence on flood extent and flood patterns. In this study, we map the flood extent, flood patterns, and geomorphological responses resulting from an extreme flood event using high resolution Sentinel-2 imagery in the low gradient, non-vegetated reaches of the Río Colorado approaching its terminus on the margins of the world's largest playa, Salar de Uyuni, Bolivia. Along with multiresolution segmentation, we apply MNDWI to separate water bodies and land, and then the ICA method is used to separate wet land and dry land for the land section derived from MNDWI. The ICA is applied to the combined imagery derived from pre-and post-peak flood data, although the post-peak flood satellite imagery is not from immediately after the flood. In addition, Sentinel-2 data are used alongside very high resolution WorldView satellite imagery to investigate the geomorphological changes resulting from the flood. The combination of MNDWI and ICA methods on Sentinel-2 MSI imagery, and the complementary use of very high resolution WorldView satellite imagery, provides fundamental information on flood characteristics and impacts, and the approach potentially could be applied in this and other settings to derive data to help calibrate and validate hydrological and hydraulic models.

Study Area
The terminus of the Río Colorado is located at the southeastern margin of the Salar de Uyuni, Bolivia, the world largest playa with an area of~10,000 km 2 ( Figure 1A-C). The Río Colorado catchment is located in the southern part of the internally-draining Altiplano Basin, which formed as part of the Andean oceanic-continental convergent margin ( Figure 1A,B). The catchment comprises upper Ordivician to Tertiary clastic sedimentary and igneous rocks, with Quaternary sediments widespread [30,31]. Despite some prominent fault escarpments in the catchment, the study area has been tectonically quiescent in the late Pleistocene and Holocene [32,33].
Bolivia, the world largest playa with an area of ~10,000 km 2 ( Figure 1A-C). The Río Colorado catchment is located in the southern part of the internally-draining Altiplano Basin, which formed as part of the Andean oceanic-continental convergent margin ( Figure 1A,B). The catchment comprises upper Ordivician to Tertiary clastic sedimentary and igneous rocks, with Quaternary sediments widespread [30,31]. Despite some prominent fault escarpments in the catchment, the study area has been tectonically quiescent in the late Pleistocene and Holocene [32,33]. The Altiplano Basin has an overall semiarid climate, with an annual precipitation of more than 800 mm in the north and less than 200 mm in the south [35] and an annual evapotranspiration potential of 1500 mm [36,37]. The north-south decrease in precipitation is due to the prevailing low pressure weather systems, whereby strong low-level northwesterly winds with warm, moist, and unstable air flow along the eastern flank of the central Andes and give rise to convective precipitation. Poleward low-level airflow helps to maintain the intense convection [38]. Owing to its southerly location, the Río Colorado is ephemeral with river flow occurring mainly in response to thunderstorms in the austral summer (December through March [39]). Although there are no flow gauging records, small to moderate (sub-bankfull) river flow events occur one or more times in most years, with larger events (bankfull or above) occurring at least once every few years. According to local accounts, following heavy rainfall and significant flow in the Río Colorado and other local rivers, water depths in Salar de Uyuni can reach up to 10 m deep, but analysis of Landsat imagery  indicates that the lake typically dries out in the intervening winter months [40]. Field data on sediment loads are limited, but grain-size analyses indicate that the lower Río Colorado system is dominated by silt and clay with subordinate very fine sand [41].

Materials
Daily precipitation data in the Uyuni area were collected for the period 1976 through 2017 from the Bolivian Servicio Nacional de Meteorología e Hidrología. Although we have precipitation data from only one gauge, the major moisture supply by strong low-level northwesterly winds along the eastern flank of the central Andes makes the precipitation pattern consistent in the drainage catchment of the study area [40].
Launched in June 2015, Sentinel-2A Level-1C (L1C) Multispectral Instrument (MSI) data are publically available for free download from Sentinels Scientific Data Hub (https://scihub.copernicus. eu/). MSI data are characterized by a span of 13 spectral bands with a spatial resolution on the ground ranging from 10 m to 60 m ( Table 1). These spectral and spatial resolutions, along with the availability of images free of charge, make Sentinel-2 very appealing for water mapping and flood monitoring as well as for other applications such as crop monitoring in agriculture (e.g., [28,42,43]). The Sentinel-2 L1C datasets are the standard product of Top of Atmosphere (TOA) reflectance and were pre-processed to generate and format bottom of atmosphere (BOA) Level-2A products using the Sen2Cor processor (version 2.3) under Anaconda Python platform [44]. Sen2Cor enables Sentinel-2 L1C products to be processed for physical atmospheric, terrain and cirrus correction and creates BOA reflectance corrected bands [43]. The output product format is a collection of JPEG-2000 images with bands reproduced at three different resolutions (10, 20 and 60 m). Sen2Cor also produces a resampled version of the higher resolution bands (10 and 20 m) to 20 and 60 m, respectively. In this study, we used bands with resolution of 20 m to extract flood extent and undertake classifications. In the SWIR bands (spatial resolution of 20 m), water absorbs solar radiation more strongly than in the visible (spatial resolution of 10 m), and these are subjected to MNDWI and ICA methods to extract flood extents in the subsequent analysis. The use of the 20 m version of the visible and NIR bands is justified by the need to avoid any artifacts that would have arisen by sharpening lower resolution SWIR bands.

Methodology
Daily precipitation data can be used for maximum daily precipitation frequency analysis and to pinpoint the peak flood events. Pre-and post-flood imagery can then be selected for extracting flood extents and patterns and analyzing the associated geomorphological changes ( Figure 2).
In this study, the daily precipitation data were processed to extract maximum daily precipitation and return periods. The maximum daily precipitation frequency analysis showed that the peak rainfall  (Figure 3), and this corresponded with an extreme flood. This flood event was within the observation periods of Sentinel-2. The Sentinel-2 image available from immediately after this event is from 9 January 2017 with cloud cover of up to 50% but the image on 19 January 2017 is less cloudy, particularly in the medial and distal regions of the lower Río Colorado, which is the main region of interest for investigating geomorphological changes. Therefore, we selected the Sentinel-2 data acquired on 19 January 2017 for the post-flood image. For the pre-flood image, the data acquired on 20 December 2016 were selected, due to high cloud cover on the image of 30 December 2016 (Table 2). This is justified by the limited rainfall between 20 December 2016 and 6 January 2017.
The selected pre-flood and post-flood Sentinel-2 data were then processed for flood extent extraction and assessments ( Figure 2). Before processing, using visual inspection, clouds and cloud shadows were detected and removed. The lower Río Colorado is characterized by low gradients (~0.000575 m/m declining to~0.000148 m/m) and a surrounding, largely homogeneous, red-brown oxidized, fine-grained surface with no vegetation, although white, salt-rich surfaces also occur locally [40,41,45]. The main region of interest-namely the medial and distal parts of the river and floodplain-was almost free of clouds, but some sporadic clouds and cloud shadows were easily identified due to their strong contrast with surrounding areas. Post-flood data were first used to separate then-existing water bodies and land according to MNDWI, along with multi-resolution segmentation thresholding. MNDWI has been reported to successfully extract water bodies because water absorption of solar radiation in SWIR bands is stronger than in NIR and visible bands [46]: where ρ green is reflectance in Band 3 and ρ SW IR is reflectance in Band 11 for Sentinel-2 [28,47]. As indicated in Section 3.1, we used Band 3 and Band 11 with 20 m spatial resolution derived from Sen2cor output files for MNDWI calculation, although Band 3 had an original spatial resolution of 10 m. The threshold is generally set to zero to extract water bodies from MNDWI, but in practice, different data acquisition conditions (e.g., satellite platforms and regions) mean that the characteristics of water index in specific images determine the threshold values. Water indices algorithms also tend to produce positively high pixel values for water bodies while other surface materials are negative, resulting in a typical bi-model of distribution of foreground and background pixels on the image histogram [28]. Histogram shape-based Otsu's method was then used to determine the threshold value for the surface water (0.24). In this study, an object-based multiresolution segmentation algorithm was first applied on MNDWI to build homogeneous polygons, and their spectral mean values were calculated [48][49][50]. Parameters of multiresolution segmentation processes were set using the values in Table 3. For multiresolution segmentation of MNDWI, only one layer (MNDWI results) was used and therefore the image layer weight was 1. We also used a thematic layer to specify regions of interest for segmentation. The scale parameter is used to determine the maximum allowed heterogeneity of the resulting image objects, while the composition of homogeneity criterion is used to define the object homogeneity to which the scale parameter refers. The polygons with the mean values greater than the threshold value were categorized as water and polygons with the mean values less than the threshold value were categorized as land ( Figure 4).       Regions of interest on land areas were subsequently created using the segmentation results of MNDWI and applied to a subset of the land area from the stacking data of pre-flood and post-flood Sentinel-2 imagery. Independent component analysis (ICA) was used on the stacked Sentinel-2 data, and flood components (flood IC) were identified through resultant independent components.
The underlying idea of independent component analysis (ICA) is that a set of multivariate signals can be decomposed into statistically independent sources without prior knowledge about the statistics of the source [25,26]. Given a set of random variables = (x1, x2, …, xn), and assuming that these variables are a linear mixture of independent components: = (s1, s2, …, sn), the model of ICA is as follows: where is the unknown matrix of mixing parameters, is an observation vector, and ICA could estimate , which is capable of computing its inverse . To extract source signals (independent components), Equation (2) is rewritten as: For a pre-and post-flood image stack, changed pixels due to moisture differences are accurately distinguished as the uncorrelated data using the algorithm ICA [24,51]. This study involved NIR and Regions of interest on land areas were subsequently created using the segmentation results of MNDWI and applied to a subset of the land area from the stacking data of pre-flood and post-flood Sentinel-2 imagery. Independent component analysis (ICA) was used on the stacked Sentinel-2 data, and flood components (flood IC) were identified through resultant independent components.
The underlying idea of independent component analysis (ICA) is that a set of multivariate signals can be decomposed into statistically independent sources without prior knowledge about the statistics of the source [25,26]. Given a set of random variables X = (x1, x2, . . . , xn), and assuming that these variables are a linear mixture of independent components: S = (s1, s2, . . . , sn), the model of ICA is as follows: where A is the unknown matrix of mixing parameters, X is an observation vector, and ICA could estimate A, which is capable of computing its inverse W. To extract source signals (independent components), Equation (2) is rewritten as: For a pre-and post-flood image stack, changed pixels due to moisture differences are accurately distinguished as the uncorrelated data using the algorithm ICA [24,51]. This study involved NIR and SWIR bands of Sentinel-2, which have a significant reduction in reflectance in the post-flood image and therefore help to identify inundated areas [52]. The pre-and post-flood Sentinel-2 imagery was then stacked to create a 12-band composite raster, which was subsequently processed using ICA transformation. Experiments on increasing complexity of the ICA, such as increasing the number of iterations and thresholds, indicated that those changes would not produce prominent changes in the resulting ICs [24]. ICA parameters were then set as in Table 4. ICA transformation created 12 independent components, which were used to further determine flooding areas. The criterion to select a flood IC is a significant reduction in reflectance of the NIR and SWIR bands due to residual soil moisture [24]. In this study, IC band 6 was selected as the flood IC for separating dry land and wet land. Using multi-resolution segmentation thresholding (Table 3) along with experiments on density slice, dry land and wet land were extracted from the flood IC. Both the segmentation of land and water as well as dry land and wet land were assessed using the visual interpretation map of the study area. Confusion matrices were used to evaluate the classification results from MNDWI and flood IC including commission error, omission error, producer's accuracy, user's accuracy, overall accuracy and Kappa coefficient.
Multi-resolution segmentation thresholding was processed using a Trimble software eCognition Developer 9.4. Sentinel-2 bands with 10 m resolution along with WorldView imagery available on Google Earth were used to visually investigate changes to channel and floodplain morphology.

Results
The main region of interest has an area of 370 km 2 (no cloud or cloud shadows included) and is located in the medial and distal parts of the lower Río Colorado ( Figure 1D). The surface is composed of red-brown oxidized, fine grained sediment and salt, and the playa is usually dry owing to strong evapotranspiration in the periods after small to moderate floods [40]. However, ephemeral water bodies could still be identified after the extreme flood event. The difference in reflectance between pre-and post-flood imagery enables water and land to be distinguished, as well as moisture saturated surfaces (wet land) and drier surfaces (dry land) using bands of NIR and SWIR.

Visual Interpretation
Based on spectral characteristics of the post-flood image, surface materials were categorized into four types: water, dry land, wet land, salty surfaces. Water surfaces are characterized by strong absorption (to zero reflectance) on SWIR bands (Bands 11 and 12), although reflectance on RGB and NIR bands is different in various locations ( Figure 5A). The edge of the playa, with reflectance typified as strong absorption at the NIR band and minor absorption in blue and green bands, indicates relatively deep and pure water. Some locations near the playa, however, have strong reflectance in the green and blue bands, indicating shallow water with salt in the bed ( Figure 5A, green line). Water pools characterized by strong absorptions on green and blue bands and elongated between channels were also categorized into this class. Dry land is characterized by strong reflectance in red and NIR and SWIR bands, but blue and green bands show absorption ( Figure 5B). By contrast, the spectral patterns of wet land are more consistent and are characterized by strong absorption both in green and blue as well as SWIR bands, whereas red and NIR bands have strong reflectance ( Figure 5C). Salty surfaces are characterized by high reflectance on red and NIR bands with minor absorption on blue and green bands but strong absorption in SWIR bands, of which band 12 is low to zero reflectance ( Figure 5D). Salty surfaces were initiated from wet land, which tends to have high water content. Also wet land and salty surfaces have the same pattern on absorption of SWIR bands. To simplify the classes and focus on mapping wet land, wet land and salty surfaces were therefore categorized into one class as wet land. The visual interpretation map based on spectral characteristics was created and subsequently used as a reference map for accuracy assessment.  Figure 1D indicate the locations for the different lines shown on each graph.

Accuracy Assessment
Using Otsu's method, MNDWI-based segmentation separated water with an area of 93 km 2 and land with an area of 277 km 2 ( Figure 6A). The overall accuracy for water separation is approximately 96% with a Kappa coefficient of 0.9 ( Table 5). The producer's accuracy and user's accuracy for both water and land were greater than 88% ( Table 5). The main errors occurred in the water-land transition zone at the edge of the playa.
Segmentation of the flood IC along with density slice experiments separated wet land using a threshold value of 0.39, resulting in wet land with an area of 142 km 2 and dry land with an area of 135 km 2 ( Figure 6B). The overall accuracy for wet land and dry land separation is 83% and its Kappa coefficient is 0.65 ( Table 6). The producer's accuracy and user's accuracy for both wet land and dry land were between 78.6% and 86.5% ( Table 6). The main errors for wet land and dry land separation also occurred in the water-land transition zone at the edge of the playa as well as in dry land-wet land transitions.

Accuracy Assessment
Using Otsu's method, MNDWI-based segmentation separated water with an area of 93 km 2 and land with an area of 277 km 2 ( Figure 6A). The overall accuracy for water separation is approximately 96% with a Kappa coefficient of 0.9 ( Table 5). The producer's accuracy and user's accuracy for both water and land were greater than 88% ( Table 5). The main errors occurred in the water-land transition zone at the edge of the playa.
Segmentation of the flood IC along with density slice experiments separated wet land using a threshold value of 0.39, resulting in wet land with an area of 142 km 2 and dry land with an area of 135 km 2 ( Figure 6B). The overall accuracy for wet land and dry land separation is 83% and its Kappa coefficient is 0.65 ( Table 6). The producer's accuracy and user's accuracy for both wet land and dry land were between 78.6% and 86.5% ( Table 6). The main errors for wet land and dry land separation also occurred in the water-land transition zone at the edge of the playa as well as in dry land-wet land transitions.

Flood Extent, Flood Patterns, and Geomorphological Changes
MNDWI-based water body mapping indicated major flow paths extending from the middle part of the study area to the playa, particularly near the terminus of the trunk channel as well as in low topographic regions between channels ( Figures 1D and 6A). The prominent water body on the northeastern corner of the study area is attributed to external flow sources separate from the Río Colorado ( Figure 6A). The water bodies in the land region with no connection to the main water body in the playa are topographic lows between abandoned channels ( Figure 6A).
Wet land and dry land distributions were consistent with major flow paths identified from the water body mapping ( Figure 6B). Larger areas of dry land were located in the proximal part of the study area, whereas wet land increased in area in the more medial and distal parts, being distributed along the major flow paths. Substantial overbank flooding occurred along the trunk channel, in part as a result of flow depth exceeding channel bankfull level, but also as a result of flow conveyance through the numerous crevasse channels and splays that characterize the lower Río Colorado [45]. The overbank floodwater flowed farther downslope, forming relatively confined flow in topographic lows between channels. In addition, wet land mapping revealed that overbank flow adjacent to the trunk channel crossed abandoned channels towards the floodplain margin (southwestern part of the dry land-wet land map), which is indicative of the extent of the overbank flow that occurred during the flood peak. Wet land mapping also indicated that a floodplain channel was reactivated with flow  Overall Accuracy 96% Kappa Coefficient 0.9 Table 6. Statistics of classification accuracy between dry land and wet land.

Flood Extent, Flood Patterns, and Geomorphological Changes
MNDWI-based water body mapping indicated major flow paths extending from the middle part of the study area to the playa, particularly near the terminus of the trunk channel as well as in low topographic regions between channels (Figures 1D and 6A). The prominent water body on the northeastern corner of the study area is attributed to external flow sources separate from the Río Colorado ( Figure 6A). The water bodies in the land region with no connection to the main water body in the playa are topographic lows between abandoned channels ( Figure 6A).
Wet land and dry land distributions were consistent with major flow paths identified from the water body mapping ( Figure 6B). Larger areas of dry land were located in the proximal part of the study area, whereas wet land increased in area in the more medial and distal parts, being distributed along the major flow paths. Substantial overbank flooding occurred along the trunk channel, in part as a result of flow depth exceeding channel bankfull level, but also as a result of flow conveyance through the numerous crevasse channels and splays that characterize the lower Río Colorado [45]. The overbank floodwater flowed farther downslope, forming relatively confined flow in topographic lows between channels. In addition, wet land mapping revealed that overbank flow adjacent to the trunk channel crossed abandoned channels towards the floodplain margin (southwestern part of the dry land-wet land map), which is indicative of the extent of the overbank flow that occurred during the flood peak. Wet land mapping also indicated that a floodplain channel was reactivated with flow draining to a topographic low between this channel and an adjacent abandoned channel, and consequently formed a major flow path to the inland water bodies. This significant overbank flooding event led to prominent geomorphological changes including crevasse channel and splay development and channel avulsion, as well as other changes to floodplain morphology such as erosion cell development.
For instance, along many reaches of the trunk channel, existing crevasse channels and splays were enlarged while levee breaching resulted in formation of new crevasse channel and splays ( Figures 7A,B and 8). Some floodplain channels also experienced crevasse channel and splay development ( Figure 7C-F).
Remote Sens. 2017, 9, x FOR PEER REVIEW 12 of 18 draining to a topographic low between this channel and an adjacent abandoned channel, and consequently formed a major flow path to the inland water bodies. This significant overbank flooding event led to prominent geomorphological changes including crevasse channel and splay development and channel avulsion, as well as other changes to floodplain morphology such as erosion cell development.
For instance, along many reaches of the trunk channel, existing crevasse channels and splays were enlarged while levee breaching resulted in formation of new crevasse channel and splays ( Figures 7A,B and 8). Some floodplain channels also experienced crevasse channel and splay development ( Figure 7C-F).  Near the terminus of the trunk channel, splay development was associated with channel avulsion (Figure 9). Avulsion involves the shift of a channel to a new position on the floodplain, and commonly involves the formation of a new channel and abandonment of the old channel. Satellite imagery shows that the newly-formed channel was initiated as separate splay channels ( Figure 9A), two of which subsequently developed an anabranching pattern before rejoining into a single channel ( Figure 9B). These channels directed flow away from the original trunk channel course, which was abandoned ( Figure 9B). During the 6 January 2017 extreme flood, one of these newly-formed channels was also abandoned to leave one dominant channel ( Figure 9C), and the avulsion is now complete.
In other locations, overbank water that was conveyed through splays promoted the development of erosion cells (Figure 10). Many erosion cells form where abandoned channels are oriented at right angles to overbank flow paths, typically leading to narrow and deep channelized sections in the scour-transport-fill sequence ( Figure 10A,B). In topographic lows between channels, channelized sections tend to be wider and shallower ( Figure 10C,D). These wider and shallower channelized sections commonly form anabranching networks, with flow paths likely having changed rapidly during the flood ( Figure 10D).  Near the terminus of the trunk channel, splay development was associated with channel avulsion (Figure 9). Avulsion involves the shift of a channel to a new position on the floodplain, and commonly involves the formation of a new channel and abandonment of the old channel. Satellite imagery shows that the newly-formed channel was initiated as separate splay channels ( Figure 9A), two of which subsequently developed an anabranching pattern before rejoining into a single channel ( Figure 9B). These channels directed flow away from the original trunk channel course, which was abandoned ( Figure 9B). During the 6 January 2017 extreme flood, one of these newly-formed channels was also abandoned to leave one dominant channel ( Figure 9C), and the avulsion is now complete. Near the terminus of the trunk channel, splay development was associated with channel avulsion (Figure 9). Avulsion involves the shift of a channel to a new position on the floodplain, and commonly involves the formation of a new channel and abandonment of the old channel. Satellite imagery shows that the newly-formed channel was initiated as separate splay channels ( Figure 9A), two of which subsequently developed an anabranching pattern before rejoining into a single channel ( Figure 9B). These channels directed flow away from the original trunk channel course, which was abandoned ( Figure 9B). During the 6 January 2017 extreme flood, one of these newly-formed channels was also abandoned to leave one dominant channel ( Figure 9C), and the avulsion is now complete.
In other locations, overbank water that was conveyed through splays promoted the development of erosion cells ( Figure 10). Many erosion cells form where abandoned channels are oriented at right angles to overbank flow paths, typically leading to narrow and deep channelized sections in the scour-transport-fill sequence ( Figure 10A,B). In topographic lows between channels, channelized sections tend to be wider and shallower ( Figure 10C,D). These wider and shallower channelized sections commonly form anabranching networks, with flow paths likely having changed rapidly during the flood ( Figure 10D).  In other locations, overbank water that was conveyed through splays promoted the development of erosion cells ( Figure 10). Many erosion cells form where abandoned channels are oriented at right angles to overbank flow paths, typically leading to narrow and deep channelized sections in the scour-transport-fill sequence ( Figure 10A,B). In topographic lows between channels, channelized sections tend to be wider and shallower ( Figure 10C,D). These wider and shallower channelized sections commonly form anabranching networks, with flow paths likely having changed rapidly during the flood ( Figure 10D).

Strengths and Limitations
This study combined MNDWI and ICA methods to extract flood extent near the terminus of the Río Colorado. MNDWI has been widely used to map water bodies, particularly for Sentinel-2 MSI data since they have SWIR bands with a spatial resolution of 20 m [28,29,53]. ICA enables identification of change classes for the stacking layer of pre-and post-flood imagery [24]. In addition, many large or extreme flood events in drylands can be attributed to intense convective thunderstorms that are characterized by spatial and temporal randomness and short durations. These hydroclimatic characteristics of drylands make field data acquisition difficult and expensive. The ICA method with no priori for classification is useful for dryland flood mapping. The ICA results do not require any post-processing such as smoothing flooding maps using clumping procedures and moving windows, which are usually used for large river systems [24]. The ICA method also resulted in precise accuracy when delineating narrow channels or dry land ( Figure 6B). The optimized combination of MNDWI and ICA thus has proven particularly useful for mapping flood extent and patterns in this nonvegetated terminal dryland river system but might also be applied to flood mapping in other hydroclimatic settings.
Owing to cloud cover, the post-flood imagery was selected from 19 January 2017, 13 days after the flood event of 6 January 2017. During the period between the flood peak and the time of postflood image acquisition, it was mostly rainy, although there were no flood events and only two separate days had daily precipitation >5 mm (19 mm and 12 mm). Hence, although the post-flood imagery was not acquired immediately following the flood peak, subsequent rain events meant that water bodies or wet land did not rapidly dry out. This study indicates the potential extension (or wider periods) for selecting post-flood imagery.

Strengths and Limitations
This study combined MNDWI and ICA methods to extract flood extent near the terminus of the Río Colorado. MNDWI has been widely used to map water bodies, particularly for Sentinel-2 MSI data since they have SWIR bands with a spatial resolution of 20 m [28,29,53]. ICA enables identification of change classes for the stacking layer of pre-and post-flood imagery [24]. In addition, many large or extreme flood events in drylands can be attributed to intense convective thunderstorms that are characterized by spatial and temporal randomness and short durations. These hydroclimatic characteristics of drylands make field data acquisition difficult and expensive. The ICA method with no priori for classification is useful for dryland flood mapping. The ICA results do not require any post-processing such as smoothing flooding maps using clumping procedures and moving windows, which are usually used for large river systems [24]. The ICA method also resulted in precise accuracy when delineating narrow channels or dry land ( Figure 6B). The optimized combination of MNDWI and ICA thus has proven particularly useful for mapping flood extent and patterns in this non-vegetated terminal dryland river system but might also be applied to flood mapping in other hydroclimatic settings.
Owing to cloud cover, the post-flood imagery was selected from 19 January 2017, 13 days after the flood event of 6 January 2017. During the period between the flood peak and the time of post-flood image acquisition, it was mostly rainy, although there were no flood events and only two separate days had daily precipitation >5 mm (19 mm and 12 mm). Hence, although the post-flood imagery was not acquired immediately following the flood peak, subsequent rain events meant that water bodies or wet land did not rapidly dry out. This study indicates the potential extension (or wider periods) for selecting post-flood imagery.
Unlike the previously-studied Colorado Front Range in the United States [24], the non-vegetated nature of the lower Río Colorado study area provides excellent conditions for satellite data acquisition and geomorphic change detection. Compared to Landsat imagery, the resolution of Sentinel-2 MSI imagery is as low as 10 m for visible to NIR bands, which improves visualization of geomorphological changes (Figures 8 and 9). Like Landsat imagery, however, the accuracy of boundary delineation of wet land and dry land for visual assessment is undermined by resolution even as low as 10 m. In addition, salty surfaces are generally characterized by strong reflectance after long dry seasons [40], but in the rainy season, strong variations from more proximal to more distal reaches occur due to different magnitudes of overbank flow. During the recession period following the flood peak, overbank flooding stops first in the proximal reaches while more distal reaches may still experience overbank flow due to the marked downstream reduction in channel cross-sectional area [45,54] and associated flow displacement. Hence, overbank water tends to dry out unevenly, resulting in salt precipitation in the more proximal reaches while reaches farther downstream may still be wet. These variations make it difficult to distinguish wet and dry salt and lead to omission errors. Errors also occur at the playa margin due to the presence of mixed pixels ( Figure 6B).

Morphodynamics Resulting from Extreme Flood Events in Terminal Dryland River Systems
The results of this study support previous studies that investigated channel-floodplain morphodynamics in the study area over longer (multi-annual and multi-decadal) timescales [39]. In particular, the results show how widespread, pronounced and rapid channel-floodplain changes occurred even in response to an individual extreme flood. During the 6 January 2017 extreme flood event, adjustments to flow sources and intensity, combined with the lack of vegetation on the fine-grained (predominantly silt, clay) sediments, led to development of crevasse channels and splays, channel avulsion and erosion cell development (Figures 7-10). The non-vegetated nature of the lower Río Colorado contrasts with many other dryland river systems where crevasse channels and splays, avulsions, and erosion cells and analogous features have been described, including floodplains and floodouts in central and south-eastern Australia [13][14][15]55,56] and seasonal wetlands in South Africa [16]. In many of these Australian and South African systems, the role of vegetation in the formation and development of these topographic features has been strongly emphasized; vegetation growth (e.g., grasses, shrubs, trees) increases hydraulic roughness and traps sediment, while the binding action of roots helps to prevent or slow erosion that might otherwise lead to levee breaching, crevasse and splay channel incision, and erosion cell initiation [15,16,55]. As a consequence, topographic features typically develop much more slowly over many years, decades or centuries, and the flow and sediment controls are harder to disentangle from the vegetative influences. Ongoing monitoring of the rapidly developing lower Río Colorado thus may provide additional unparalleled opportunities to document in greater detail the erosional and depositional dynamics involved in the development of such topographic features, and how these features link into wider cascades of flood-driven, channel-floodplain changes (Li et al. in review [57]).

Conclusions
The combination of MNDWI and ICA methods on Sentinel-2 MSI imagery has led to advances in flood mapping along the non-vegetated reaches of the lower Río Colorado approaching its terminus on the margins of Salar de Uyuni, Bolivia. Maximum daily precipitation frequency analysis based on a 42-year record of daily precipitation data (1976 through 2017) indicated that an approximately 40-year event (40.7 mm) happened on 6 January 2017 in the study area, and this was associated with an extreme flood. Sentinel-2 data acquired after this extreme flood were used to separate water bodies and land using modified normalized difference water index (MNDWI), and then independent component analysis (ICA) was applied to the land section of the combined pre-and post-flood images to extract flooding areas. The Río Colorado dryland river terminus system was classified into three categories: water bodies, wet land and dry land. The results were in agreement with visual assessment with an overall accuracy of 96% and Kappa of 0.9 for water-land classification and an overall accuracy of 83% and Kappa of 0.65 for dry land-wet land classification. The flooding extent mapping revealed preferential overbank flow paths, which were closely related to geomorphic changes. Sentinel-2 images along with WorldView imagery from Google Earth revealed that widespread overbank flooding resulted in enlargements of existing crevasse channels and splays and the formation of new splays. An avulsion near the terminus of the trunk channel was completed during this extreme flood, while in other locations, overbank flow was associated with erosion cell development. In this setting, adjustments to flow sources and intensity combined with the lack of vegetation on the fine-grained (predominantly silt, clay) sediments leads to unusually rapid erosion cell development, providing new insight into the development of these floodplain topographic features.