A method to detect abrupt shifts in river channel position using a Landsat‐derived water occurrence record

The lateral migration of river channels is an important control on the evolution of alluvial fans, deltas, and floodplains. Lateral migration consists of both gradual riverbank migration and abrupt shifts in location due to avulsions or cutoffs. Methods exist to measure bank migration, but abrupt shifts in position are rarely considered or are not emphasized. Here we describe a new method using Landsat‐derived water occurrence images that primarily focuses on detecting when a channel has abruptly shifted position, either from avulsion or cutoff. The method does not assume any a priori model of channel geometry or evolution. Within a given area of interest, binary channel images created from the fluvial water occurrence record are stacked through time. Then a channel shift intensity, Ai , is created by estimating the number of possible ending times for fluvial water voxels (a point in three‐dimensional space) in the stacked occurrence record. The number of possible end‐times for fluvial water voxels within a given region of the occurrence record reveals the likelihood that a reach of a river underwent an abrupt channel shift during the observation period. We present the results of this analysis for a 194 481 km2 region of the Amazonian basin. Follow‐up validation found three avulsions and 270 cutoffs within regions identified by this method. We show that areas above a threshold Ai contain an avulsion or cutoff with high probability. This highlights the method's potential to detect and quantify abrupt channel shifts at the basin scale. The method also successfully distinguishes between abrupt channel movement and complex braiding behaviour. As this method is applicable to any binary water‐masked time series images, future applications of this method have the potential to provide insight into the controls and spatial variance of avulsions and cutoffs across a variety of scales.

exist to measure bank migration, but abrupt shifts in position are rarely considered or are not emphasized. Here we describe a new method using Landsat-derived water occurrence images that primarily focuses on detecting when a channel has abruptly shifted position, either from avulsion or cutoff. The method does not assume any a priori model of channel geometry or evolution. Within a given area of interest, binary channel images created from the fluvial water occurrence record are stacked through time. Then a channel shift intensity, A i , is created by estimating the number of possible ending times for fluvial water voxels (a point in three-dimensional space) in the stacked occurrence record. The number of possible end-times for fluvial water voxels within a given region of the occurrence record reveals the likelihood that a reach of a river underwent an abrupt channel shift during the observation period. We present the results of this analysis for a 194 481 km 2 region of the Amazonian basin. Followup validation found three avulsions and 270 cutoffs within regions identified by this method. We show that areas above a threshold A i contain an avulsion or cutoff with high probability. This highlights the method's potential to detect and quantify abrupt channel shifts at the basin scale. The method also successfully distinguishes between abrupt channel movement and complex braiding behaviour. As this method is applicable to any binary water-masked time series images, future applications of this method have the potential to provide insight into the controls and spatial variance of avulsions and cutoffs across a variety of scales. Migrating alluvial river channels influence sediment erosion and deposition, and therefore affect alluvial stratigraphic architecture and floodplain evolution and morphology (Leopold, 1994;Mackey & Bridge, 1995). Channel migration can be continuous or abrupt.
Continuous migration occurs when rivers meander and the length of the centreline grows gradually (e.g. Sylvester et al., 2019). This meandering process creates floodplain on the inner point bar and erodes the floodplain on the outer cut bank. This continuous migration helps shape the landscape, while also providing vital bio-geochemical exchanges necessary to sustain complex riparian ecosystems (Latterell et al., 2006) and influencing carbon excavation and burial along riverbanks in the rapidly thawing Arctic (Douglas et al., 2022). But continuous migration is interrupted by abrupt shifts when rivers avulse to new positions on the adjacent floodplain, or when rivers migrate into themselves, causing a cutoff. These abrupt shifts are also crucial for building floodplains, because avulsions create significant sediment deposition during the process of channel establishment (Lombardo, 2017;Valenza et al., 2020), while cutoffs create significant ecological refugia and space for sediment accumulation (Constantine et al., 2014). Avulsions also create significant hazards because, during this process, there are often devastating floods (Slingerland & Smith, 2004); in 2010, avulsion-related flooding along the Indus River displaced tens of millions of people (Syvitski & Brakenridge, 2013). Remote-sensing methods that can differentiate between continuous and abrupt migration are necessary to build databases that will allow the community to investigate how and when avulsions and cutoffs occur. These databases would improve our understanding of rivers, floodplains, and associated hazards.
Progress has been made in using the remote-sensing record to measure continuous lateral channel migration, with the emergence of modern remote-sensing toolboxes focusing on tracking channel centrelines or banks (Jarriel et al., 2019;Legg et al., 2014;Monegaglia et al., 2018;Rowland et al., 2016;Schwenk et al., 2017). However, current methods that track centrelines or banks have not demonstrated strong efficacy in quantifying abrupt shifts caused by avulsions and cutoffs, mostly because existing methods focus on tracking smoothly migrating features instead of the discontinuous jumps that characterize avulsions and many cutoffs. Current solutions use single images or differences between a pair of times (Monegaglia et al., 2018;Rowland et al., 2016;Schwenk et al., 2017), rely on heuristic thresholds or a priori models of river behaviour (Monegaglia et al., 2018;Schwenk et al., 2017), or apply only to simple river morphologies (Fisher et al., 2013;Legg et al., 2014).
Developing methods to detect avulsions and cutoffs across basinto continental-scale areas would add to our understanding of fluvial processes. Regional-scale datasets of manually compiled, remotely sensed avulsion observations have already led to insight on path selection (Edmonds et al., 2016 and crevasse splay production on avulsing rivers (Lombardo, 2017). Additionally, these smaller datasets suggest that avulsive behaviour is related to stream morphology and systematically changes downstream (Valenza et al., 2020), possibly in response to downstream fining (Dingle et al., 2020). Expanding these efforts to the continental scale using the remote-sensing record could help provide answers to fundamental questions concerning the nature of both avulsions and cutoffs. To accomplish this requires reliable tools that can extract data from the satellite record in an unbiased way. Large, continental-scale datasets of field observations of avulsions and cutoffs would allow us to observe the relationship between morphology, vegetation, lithology, and anthropogenic activity on abrupt channel-shifting events. To detect avulsions and cutoffs from the remote-sensing record without significant human verification effort, a technique needs to (1) auto-adapt its threshold for abrupt channel shifts to local migration conditions and morphologies, (2) be able to process imperfect remote-sensing data, and (3) provide an unbiased indicator that a channel has abruptly shifted position in a way that deviates from the continuous migration patterns of that section of the river. A method that satisfies these three tasks is critical for automating the detection of abrupt channel shifts over large spatial areas containing many rivers of multiple scales and morphologies.
In this work, we present a method that uses the full temporal record available in a remote-sensing dataset. Rather than tracking the change in conventional channel attributes over time, the method assumes no a priori model of channel behaviour. The method can do this by building on the previous work of Pekel et al. (2016) to identify which bodies of surface water are channelized. For each body of fluvial surface water, we create a spatial-temporal stack of images where individual pixels are tracked upwards through time (using an auto-adapting threshold) to detect when a portion of a river has likely changed position. Searching through spatial-temporal space, rather than compressing the full record down to a two-dimensional (2D) occurrence frequency image, is necessary to distinguish rapid continuous migration from abrupt channel shifts. Further, we demonstrate that introducing temporal information also de-emphasizes small artifacts in spatial information, allowing us to incorporate noisy, imperfect data. While we use remote-sensing images as our input, this method is applicable to any binary masked time series of images, including laboratory experiments. This approach can be viewed as complementary to other methods that provide finer-grained detail about river behaviour but do not necessarily provide information about avulsive activity or cutoffs. The following text describes the method in full detail and presents the results of the method applied to a 441 Â 441 km region in the Amazonian basin. We demonstrate that this method can distinguish between abrupt channel movement and gradual, complex braiding behaviour. Given that the scale of the initial results is on par with recently published analyses of river dynamics (Monegaglia et al., 2018;Schwenk et al., 2017), the potential for the technique to be applied at larger spatial scales to generate unique databases of abrupt channel-shifting events is promising.

| LITERATURE REVIEW
The problem of detecting avulsions and cutoffs seems like a simple extension of existing methods that use remote sensing to study river migration by tracking centreline deformation (Jarriel et al., 2019;Legg et al., 2014;Micheli et al., 2004 Another style of analysis tracks areal changes in a river on a pixelby-pixel basis (Peixoto et al., 2009;Rowland et al., 2016). Recently, the SCREAM channel detector (Rowland et al., 2016) used a pixel-bypixel method to outline an approach for avulsion and cutoff detection using patterns of erosion and accretion obtained from the difference between binary masks. However, this approach has not been demonstrated in the literature beyond a single example (Rowland et al., 2016). Because the SCREAM approach to avulsion and cutoff detection is also threshold-dependent, it cannot be automatically applied to the diversity of river morphologies and behaviours without a systematic way to adapt that threshold. Finally, there is qualitative evidence that the channelized response variance (CRV) recently published by Jarriel et al. (2019) can detect avulsions. However, using the CRV to isolate the signature of avulsion events from other types of change has not been attempted.
Most of the existing methods to quantify channel migration rely on measuring attributes of the river obtained from an image or image mask taken during a single time step, or differencing images from two time steps. In other words, the base objects for analysis (such as centrelines or pixel-by-pixel changes) are derived from static images or differenced images. The CRV is a notable exception (Jarriel et al., 2019) that measures how channelization evolves over time and is an integrated measure over a set window of observations at a location. In this sense, the CRV incorporates more time information than most other metrics for river planform change, and is an example of how incorporating temporal information yields additional information about channel dynamics.
Preserving the temporal information of each layer in a stack of images may be more advantageous for studying both continuous channel migration and abrupt shifts away from typical migration patterns. For instance, in the case of detecting the presence of an abrupt shift in channel position, the details of the spatial change in channel position are easier to interpret when there has been a marked change in the temporal patterns of migration. Creating a temporal stack effectively de-emphasizes local noise or artifacts in spatial information, and allows us to incorporate imperfect data and metrics that track time-integrated behaviour to detect abrupt changes.

| Image acquisition, binarization, and stacking
To begin the avulsion and cutoff detection process, our method requires an ordered time series of binary water masks. In this study, we present an application of the water occurrence record recently created by Pekel et al. (2016). This is a publicly available dataset derived from the Landsat record from the years 1984 to 2018. The dataset has water occurrence available monthly and yearly, with 30 m pixels classified as no data, not water, seasonal water, or permanent water. We obtained the yearly occurrence records from Google Earth Engine and binarized them by casting the 'no data' and 'not water' pixels to 0 and the 'seasonal water' and 'permanent water' to 1 (e.g. see Figure 1).
One of the primary challenges in using the occurrence dataset is dealing with apparently disconnected stream pixels that do actually represent connected streams. Pixels may disconnect because of data quality, or because low-water conditions on small streams narrow the channel beneath the detection limit of the Landsat data and disconnect the pixels. For example, it is easy for a human to visually conclude that the eastern branch in Figure 1 still exists in panel B, since it appears in earlier and later years (panels A and C). An apparent or temporary loss of spatial connectivity can be resolved by connecting the images through the time dimension in addition to the spatial dimensions.
Our initial approach to overcome this issue was to use flattened image stacks as the base object of analysis to determine if an avulsion or cutoff occurred based on occupancy differences. Flattening a temporal image stack of water occurrence images compresses all temporal information into a single value, creating a composite image (see Figure 8c for an example of this). This technique works for certain cases, but it is problematic because the continuous lateral channel migration in many channels smears the main channel body, creating similar occupancy values for the avulsing channel and the main channel. Another issue with flattened images is that channels frequently avulse or cut off by reoccupying previous channel positions (e.g. Valenza et al., 2020Valenza et al., , 2022, or avulse and then laterally migrate over the old channel position. In such cases, flattened occupancy maps cannot distinguish between these processes. Instead, we preserved the temporal information in our stacked water occupation images, connecting them through time by creating three-dimensional (3D), 26-neighbour connected components (i.e. a 3 Â 3 Â 3 voxel cube). Both the space (x, y) and time (z) dimensions are used to connect fluvial surface water pixels in these stacked images. In this view, each pixel of an occurrence image becomes a voxel within the component (see example in Figure 2). This simple technique allows us to treat the 3D, connected components of timestacked binarized images as the object of interest while bypassing the segmentation difficulties created by noisy years in the Landsat record.
Note that while the stacking technique greatly increases the robustness of associating voxels from years with intermittent data to fluvial F I G U R E 1 An example of the image binarization performed on the Pekel et al. (2016) dataset for (a) 1996, (b) 1997, and (c) 1998 on an 81 km 2 Landsat image. White and black pixels are classified as water and not water, respectively. The portion of the stream pictured in this image is located at approximately (À17.26 , À64.51 ) and splits into two branches. Discharge variability results in apparent changes in connectivity of the pixels, especially on the east branch.
bodies, it still fails in the case where there is total data loss for a given year, or the case where there are multiple years of extremely sparse data.
Another advantage of creating 3D images with time as an axis is that it highlights the essentially ribbon-like nature of many streams as they progress through time (Figure 2). An avulsion or cutoff can be thought of as a tear or split in this ribbon-like structure, as the channel abruptly changes position. In Figure 2 there are avulsions that are clear as thinner streamers that 'tear away' from the ribbon formed by the main channel.
To demonstrate this method, we created connected components for every body of surface water within the test area used for the study. This test area was a 441 Â 441 km region in Bolivia ( Figure 3).

| Selecting fluvial surface water
Because our test input data do not distinguish between fluvial and non-fluvial surface water, we must filter out the non-fluvial components (e.g. reservoirs or lakes). To do this, we create a 2D binary mask showing the position of all fluvial water bodies using the deep learning-based fluvial classifier 'deepriver' as applied to 1 year (1998) of Landsat imagery (Isikdogan et al., 2018). This produces an image where each pixel has a value representing the likelihood that it is near the centreline of a fluvial body. This image is thresholded to remove most of the non-fluvial surface water. A final, manual masking step was applied at a coarse scale to remove remaining non-fluvial pixels over the 194 481 km 2 study area in question.
The resulting 2D binary mask containing all fluvial water pixels is then compared to the 3D components generated above. A 3D component is considered fluvial if it intersects at least one fluvial water pixel in the 2D binary mask. This filters out non-fluvial surface water as well as smaller streams that would be unreliably represented in our test input data, while leaving enough streams available for analysis.
F I G U R E 2 An example of a time-stacked channel created by taking the 3D, 26-neighbour, connected components for fluvial surface water obtained from Pekel et al. (2016). This figure is a temporally stacked record that encompasses the occurrence data shown in Figure 1. Avulsions that have occurred along this section of the channel are highlighted by arrows.

F I G U R E 3
The study area used to test the method: a 441 Â 441 km region in the Amazonian basin in Bolivia. Previously observed avulsions (Edmonds et al., 2016) are marked by the crossed circles. Avulsions found by the study method are marked with green circles; avulsions that were masked by the fluvial surface water selection protocol (see Section 3.2) are highlighted by purple circles. Light blue lines are Global River widths from the Landsat database (Allen & Pavelsky, 2018a, 2018b. WGS 84 coordinates are overlain on base map imagery.
[Color figure can be viewed at wileyonlinelibrary. com]

| Measuring the channel shift intensity A i
After masking out non-fluvial components, we measure a channel shift intensity, A i , for arbitrarily sized (we use 81 km 2 ) regions of interest (ROIs) within the study area ( Figure 4). Conceptually, A i measures the number of unique temporal terminations within a connected reach of a water-body volume. Multiple, unique temporal terminations (called end-times) can result from abrupt relocations due to avulsions or cutoffs. We measure this by tracing upward-directed pathways through the connected volume of fluvial water voxels and finding the end-times of each path. Measuring A i is a three-step process.
First, for a given starting voxel (called V start ), we find all possible end-times originating from that location ( Figure 5). V start is defined as the earliest water occurrence voxel for a given planform pixel. From V start , we trace all possible upward-connected paths, using a temporal step size of 1 year and a spatial step size of l. Visually, if the initial volume was a cube, these paths would collectively trace out a pyramid with the point originating at V start . Tracing continues upward through each transition voxel (V trans ), defined as a voxel that has at least one voxel above within step size l ( Figure 5). Once there are no more possible upward voxels, the path ends and the last voxel is defined as an end-voxel (V end ); V end is a voxel with no neighbouring voxels above it in time for the current l. Voxels that sit at the edge of channel volumes are not end-voxels unless there are no voxels that sit above them in time. Because the volume is a time stack, the z-values of each V end record the year that path terminated (called the end-time).
Finally, all unique end-times are then assigned to that V start .
Second, this process is then repeated for all possible V start in the largest connected fluvial component in the ROI (Figure 4). The A i for the ROI using a given l value is recorded as the number of unique end-times among all V start . For example, a simple channel that migrates laterally has one unique end-time (Figure 5a), whereas an avulsion creates two unique end-times (Figure 5b). Note that, regardless of the number of V start , if all have the same end-times, the sum of unique end-times for the ROI will still be the same as for any individual V start .
Third, we find all unique end-times within an ROI by repeating steps 1 and 2 and progressively increasing the length of l by 1.74 ( ffiffiffi 3 p ) until A i for the ROI is stable for five consecutive values of l (Figures 4   and 6). This stable A i is recorded as the channel shift intensity for the ROI (Figure 6). This iterative process of progressively increasing l to find a stable A i is the key feature allowing our method to adapt to variable rates of continuous channel migration, different channel widths, and to distinguish most channel braiding patterns from abrupt channel-shifting events.
When no abrupt channel shifts have occurred and the data are not noisy, A i will be 1 (e.g. see Figure 5a). Consider another idealized channel cross-section where A i is >1 due to an avulsion that occurs at year 6 ( Figure 5b). In the cross-section, it is apparent that several starting points will have possible end-times at both year 6 and year 10, and that these multiple end-times will persist for a range of l. This F I G U R E 4 Graphical depiction of the method used to measure the channel shift intensity, A i . V tr and V st refer to V trans and V start , respectively.
leads to an A i of 2. Theoretically, any system with an A i of 2 or more should host an abrupt channel movement; in practice, noise in the data causes minor artifacts. A higher value of A i instead indicates a greater likelihood that a ROI contains one or more avulsions or cutoffs. For this test analysis, an 81 km 2 (9 Â 9 km) ROI is selected and shifted across the study area 0.3 km at a time. This discrete sampling approach is an approximate way to associate abrupt channel-shifting events with distinct ROIs, and a shift of 0.3 km is smaller than the smallest avulsion length ($1500 m) measured by Edmonds et al. (2016), which suggests that we will capture most events in the test F I G U R E 6 Detected number of unique end-times A i for a single ROI as a function of the step length l. Once A i has been stable for five values of l, the stable value is declared the A i for that ROI. In this example, the ROI is assigned an A i of 10. [Color figure can be viewed at wileyonlinelibrary.com] area. A i values will be greatest for ROIs that directly overlie the centre of active, temporally branching regions, and will decline laterally as the ROI moves away from these regions. For these parameters, the current analysis took $36 h to run over the 440 km 2 study area on a Google Cloud compute instance running 12 vCPUs and 32GB of memory. This implies that with enough computational resources and/or further optimization of the method, much larger areas could be analysed within reasonable timeframes.

Maps of sampled
3.5 | Using A i maps to measure channel-shifting activity for different regions within the study area

| Test cases
To demonstrate our method and test its performance, we applied it to several synthetic test regions (Figures 7a-d), as well as real water occurrence data (Figures 7e-h). The method performed as expected for the synthetic data, where a static channel and a channel with a constant migration rate have an A i of 1 (Figures 7a and b) and synthetically generated avulsions and cutoffs have A i values of 2 (Figures 7c and d). The method also performs as expected on four real-world test samples (Figures 7e-h). A braided river without avulsions outside of the braid belt receives an A i of 1 (Figure 5e). In sample two, one large avulsion off a rapidly migrating channel creates an  region lead to different values of A i . Sample three contains two avulsions and has an A i value of 3 (Figure 7g). In sample four, the channel section has avulsed and cut off multiple times and has a high A i of 7 (Figure 7h).

| Generating an A i map for a region of the Amazonian basin
After testing across the cases in Figure 7, we generated a map of A i for all ROIs within the entire study area (Figure 3). In the A i map, the channel network is interspersed with areas of low and high A i This shows that regions with high values of A i are likely to contain an abrupt channel-shifting event.

| A i threshold analysis
We then performed a more systematic examination of areas of high A i . Regions with a higher A i value contain a greater number of endtimes, and as such are more likely to contain abrupt channel movement events (Figure 9). In general, as a chosen threshold A i value (see Section 3.5) decreases, the probability of an ROI containing an event goes down (Figure 9a) and the total area of the regions above that threshold increases (Figure 9b). For the study area in question, it is clear that for A i > 7, the method will reliably identify regions with a 100% probability of having undergone an abrupt channel-shifting event (Figure 9a). Even below this threshold, the probability that a region with an A i value of 2 contained an avulsion or cutoff was >50%. The specific A i value of the 100% probability threshold found in this study may differ between study areas or datasets, and so users should begin with the highest-certainty areas and evaluate successively lower threshold values, as by definition each set of regions of a particular A i value contains the set of all regions with higher A i values ( Figure 10). The density of abrupt channel-shifting events as a function of the threshold A i value was also evaluated (Figure 9c). In general, as the threshold A i value decreases, the event density decreases (Figure 9c). For the highest A i values, the system approaches a density of two avulsions or cutoffs per 100 km 2 .

| Detected avulsions and cutoffs
To evaluate the correspondence between A i and the occurrence of avulsions and cutoffs, all detected regions for each threshold value were evaluated by human inspection for the number of avulsions and cutoffs they contained. Avulsions and cutoffs were evaluated for each region and tallied manually by two or three people to achieve a degree of objectivity. The whole study area was found to contain three avulsions and 270 cutoffs for all regions above an A i threshold of 2.
We first evaluate the method's performance with regard to avulsions. We compared the avulsions we located in our test dataset to those found in the same area by Edmonds et al. (2016). Using manual inspection of Landsat timelapse data, those authors found an additional eight avulsions beyond those identified by our method.
Videos of both the previously found avulsions and the avulsions found by the method can be found in the Zenodo dataset. While there is overlap between the previously observed avulsions and the avulsions  F I G U R E 1 0 Yellow pixels depict regions above the current threshold A i value over the entire study area. As the A i value is progressively lowered, the number of regions above the threshold increases. Each pixel is 81 km 2 . [Color figure can be viewed at wileyonlinelibrary.com] found by the method (see Zenodo dataset for a marked map of avulsion locations), our method struggles with detecting avulsions that occur on scales approaching the pixel resolution of our input data. This is because these avulsions are small enough to be masked out of the test input dataset used in this example. The method also struggles with avulsions that are large enough to surpass the 9 Â 9 km ROI over which each A i value is measured. Increasing this window size could help capture these very large events. The size of these events makes them the easiest to identify manually, however, and therefore the likeliest to have been previously identified. The method also struggles with avulsions when the surface water data is exceptionally sparse or noisy. Interestingly, the highest A i region is not an avulsion, but a ROI where numerous, large cutoffs have occurred frequently in a small area.
The method performed much better with cutoff detection, as validated by manual inspection. Above an A i of 7, 100% of regions contained cutoffs and there were only three regions with an A i of 1 where cutoffs occurred but were not detected. As such, while the method most often over-detects abrupt channel-shifting activity (as only $60% of regions with an A i of 2 contain events), there do exist rare cases where a region underwent an abrupt shift (i.e. cutoff or avulsion) that was not detected. This type of error (undetected events) occurs in $3% of cases, and appears if cutoffs occur at scales below the length scale associated with evaluating A i at a minimum number of consecutive step lengths, l, until A i reaches a stable value.
In this study, A i needed to be stable for five consecutive values of l.
This stability criterion means that, at a minimum, A i represents the number of distinct end-times evaluated at an l of $200 m with larger values of l often being necessary to ensure the stability of the A i measurement for that region. If a small cutoff happens below this scale, it is possible that it will go undetected and the A i for the region in question will be 1. Conceptually, this error is a misclassification of an abrupt change as being continuous due to too small a spatial change between measurement time steps. This error is a trade-off inherent in having a stability criterion that is dependent on the step length. In general, this dependency is desirable as it ensures that braiding activity within the channel belt is mostly ignored and allows the method to adapt automatically to variable channel widths and migration rates.

| DISCUSSION
The A i metric successfully and objectively detects abrupt channelshifting events as A i increases (Figure 9), and quantifies the spatial distribution of abrupt channel-shifting events across noisy data that contains rivers of differing scales, morphologies, and change intensities ( Figure 8). The flexibility and noise robustness of our method also mean that while higher A i values indicate higher probabilities or quantities of abrupt channel-shifting events, the magnitude of A i is not related to any measure of the size of these events and does not specifically discriminate between avulsions and cutoffs. For example, a region with a high A i could be one avulsion and several cutoffs or, occasionally, simply a noisy measurement with one avulsion. Figure 11 demonstrates the relationship between A i and morphology for typical low-, medium-, and high-A i ROIs. Low-A i regions will generally have only one or two avulsions or cutoffs. Medium-A i regions may have several events, while high-A i regions are generally highly active in terms of abrupt channel movement events, regardless of whether these events are avulsions, cutoffs, or both. Importantly, however, there is no general trend of higher-A i areas appearing more or less braided; the method correctly discriminates relatively uncommon abrupt channel movements from the more common gradual braiding or meandering observed in much of the study area.
Due to this imperfect correspondence between A i and the magnitude or nature of any individual abrupt channel-shifting event in an ROI, A i is best treated as a probability of the presence or absence of abrupt channel-shifting events within a region. Further work could attempt to increase the predictive power of A i on the count or nature of abrupt channel events, including by an advancement of the analysis method.
We compared the avulsions we located in our test dataset to those found in the same area via manual inspection by Edmonds et al. (2016). In this area, those authors found eight additional avulsions beyond the three we found. This is because those authors used human inspection of Landsat timelapse data, and we solely used the A conservative estimation of fluvial pixels leads to some abrupt channel-shifting events not being detected because they were incorrectly classified as occurring in non-fluvial surface water in the method. In spite of this, A i maps of both the unmasked surface water data and the surface water dataset masked using deepriver show similar spatial patterns for the fluvial areas that are present within both maps.
It is important to note that the method is insensitive to branching, braiding, and continuous lateral channel migration. This insensitivity to braiding and lateral migration arises because the requirement that A i be stable for a range of l creates an adaptive threshold for dividing continuous from abrupt channel activity. Examples of this can be seen in both the synthetic constantly migrating stream as well as in the braided test case (see Figures 7b aznd e). Branching streams do not have a high value for A i because they branch in space and not in time, and our method is designed to ignore spatial branching patterns.
A i is an accurate indicator of abrupt channel shifting within a region if the appropriate threshold value for A i is chosen (100% detection rate), but choosing this threshold is not straightforward. Figure 9b indicates that for future studies this threshold could be selected automatically by finding where the growth rate in the total area occupied by detected regions above a threshold A i goes past a certain value. Defining a threshold A i has the potential to be useful for creating continental-scale databases of abrupt channel-shifting events from the remote-sensing record that have previously been challenging to create. Such databases will allow the community to better examine how factors such as floodplain morphology, sediment type, and rate of continuous lateral channel migration affect how and when avulsions and cutoffs occur. Given that these are matters of active debate, we can only benefit from more observation.

| CONCLUSIONS
We have presented a method for estimating the spatial occurrence of abrupt channel-shifting events using Landsat-derived water occurrence observations spanning a 35-year time period. This method uses a newly created metric, A i , paired to a sampling protocol that allows for an automated and comprehensive detection of avulsions and cutoffs within a time series of a given area. The technique automatically scales itself to detect avulsions and cutoffs across a variety of spatial scales and river morphologies. This technique has been tested on an $200 000 km 2 portion of the Amazon basin and has shown the ability to detect abrupt channel-shifting events that occurred within the test area. Low rates of non-detection were observed, with $3% non-detection for 170 events. A thresholding technique was also demonstrated that allows a user to select regions of interest within the study area depending on threshold A i values. In the study area, 100% of regions above a threshold A i value of 7 contained abrupt channel-shifting events. The high fraction of regions containing avulsions or cutoffs at all threshold A i values >1 provides strong support that A i is an objective, quantitative method for detecting abrupt channel-shifting events without a priori knowledge or assumptions of fluvial geometry or behaviour.
Application of this method will allow for future exploration of variations in the spatial frequency of avulsions and cutoffs across basins.
By mapping the spatial occurrence of such events for a given time duration, the method can be used to look at how frequently channels shift abruptly within different areas of a channel belt, floodplain, fan, or even continent relative to one another. This capability has the potential to be useful in hazard analysis. Even more fundamentally, the method has the potential to assess how many avulsions and cutoffs occurred over the past 35 years throughout large sedimentary basins across the globe.