Salt Marsh and Tidal Flat Area Distributions Along Three Estuaries

Estuarine landscapes form through interactions between fluvio-coastal processes and ecological processes within the boundaries imposed by hard substrate layers and man-made dikes and dams. As estuaries are ecologically valuable areas, monitoring and quantification of trends in habitats is needed for objective comparison and management. However, datasets of tidal flat and saltmarsh habitats along entire estuaries are scarce. The objective was to compare trends of biogeomorphological areas and habitat transitions along three estuaries in the Netherlands and assess whether these are generally comparable or mainly determined by system-specific histories. We present data for these estuaries obtained by automated classification of false-color aerial imagery. The automated method allows objective mapping of entire estuaries at unprecedented resolution. The estuaries are dominated by subtidal areas and tidal flats. The tidal flats have similar area along the estuaries while saltmarsh area decreases. Collective lengths of ecologically important transitions between saltmarsh, low-energy tidal flats and water differ more between the estuaries. These variations are due to presence of mid-channel bars and shore-connected embayments. Saltmarsh area is mainly determined by the different formation and embankment histories of the estuaries, embayments and side-branches. Much of the past saltmarsh flanking the estuaries was lost due to past land reclamation. In one system, ecologically important low-energy tidal flats are reduced by a sudden decrease of tidal amplitude, causing increase of subtidal area at the cost of intertidal area. Large areas of high-energy tidal flats in one estuary remain unexplained. The automated method can be applied in other estuaries, provided that high-quality areal imagery is available. Extensions of the data to other estuaries would allow for system-scale trend comparison between estuaries of ecologically relevant biogeomorphological characteristics.


INTRODUCTION
Estuarine landscapes form through interactions of tidal currents, wind-generated waves and river discharge and the balance of sediments transported by these processes, all of which may be modified by biological processes. At the same time, estuaries are topographically forced by antecedent morphology, especially by constraining, hard or unerodible substrates (Townend, 2012;de Haas et al., 2018;de Haas et al., 2019), henceforth referred to as initial conditions. In addition, estuaries are also directly and indirectly, through reduction of tidal prism, constrained by the effects of land reclamation, embankments, barrages and storm-surge barriers. Constraints on width and on tidal amplitude immediately modify the hydromorphodynamics and gradually change the long-term, large-scale development of the entire system in multiple ways. The seaward boundary conditions of tidal hydrodynamics, wave climate and sediment fluxes are not only directly limited by constraints at the estuary mouth, but also by landward reduction of shallow areas, which reduces the flood storage. This, in turn, reduces the tidal prism in the entire system (Friedrichs and Aubrey, 1988) and may, on the longer term, lead to shrinking of the estuary (Braat et al., 2017(Braat et al., , 2019de Haas et al., 2018). This is also partly determined by river discharge and fluvial sediment supply (de Haas et al., 2018). There is evidence that the sedimentary filling of late-Holocene estuaries proceeds in a large-scale pattern: the main channel in shallow estuaries generally converges as a result of approximate balancing of tidal wave height and energy loss through friction (e.g., Savenije, 2015) and balancing ebbdirected and flood-directed sediment transport (e.g., Dronkers, 2017). Hereby the embayments and reaches that are wider than the converging shape are zones of bar formation, tidal flat sedimentation and, in the most sheltered areas, saltmarsh formation (Wang and Temmerman, 2013;Leuven et al., 2018a). This results in a diverse landscape with spatial variation in tidal flat area, substrate composition and biomorphodynamics, where tidal bars (or shoals) in the estuary are expected to be more dynamic than the partially sheltered, shore-connected flats and saltmarsh flanking the estuary.

Problem Definition
Due to the large spatiotemporal variability in conditions and substrates, shallow estuaries are known to be of special ecological significance. Estuaries have high primary production, spawning zones for fish, foraging areas for birds, gradual water-land transitions, and have a unique salinity gradient from a polyhaline to mesohaline to freshwater tidal environment on the fluvialtidal transition, which, in turn, is important for certain fish species and leads to a rapid change in biodiversity of freshwater species and salinity-tolerant species (e.g., Whitfield et al., 2012;Greenberg, 2013;Telesh et al., 2013). Not surprisingly, estuaries show many interactions between hydromorphological processes and ecological processes that potentially change the landscape in morphology, sedimentology and spatial structure. Particularly eco-engineering species play a major role as they change the hydrodynamics, sedimentology and morphology in their own environment and facilitate other species (Jones et al., 1994). In particular, higher intertidal and supratidal zones are often pioneered and covered by saltmarsh species in temperate regions (e.g., Wang and Temmerman, 2013;Schwarz et al., 2018) and by mangroves in subtropical and tropical regions, while lower intertidal and subtidal zones are populated by a wealth of animal benthic species.
This leads to the question to what degree shallow estuaries share a general pattern of proportions of tidal flats and salt marsh habitats across and along the system, or whether all have unique, dominant patterns that are determined mainly by the constraints of geological setting and set by past and present human interferences. Spatial heterogeneity of, and gradients between, habitats is deemed important for many species and ecosystem quality in general (Wolter et al., 2016), but it is unknown whether the heterogeneity differs between estuaries because that requires a method for mapping and quantification that is presently lacking. Furthermore, it is at present unknown whether shore-connected tidal flats and saltmarshes and those on shoals differ in area and character, while it can be hypothesized that the shore-connected areas are less dynamic. Although methods have been developed to map benthic habitats (e.g., https://data.gov. uk/dataset/1707e638-6a2d-48f5-a534-1db0b240cc37/casi-andlidar-habitat-map) and coastal land cover (e.g., Enwright et al., 2020), datasets are insufficiently available as is evident from the fact that review papers have had to focus on local case studies (e.g., de Jonge, 2000;Whitfield et al., 2012;Telesh et al., 2013;Dashtgard and Croix, 2015) and stress the need for combined terrestrial-aquatic ecosystem studies with remote sensing techniques (Soininen et al., 2015). We only found two studies for the Netherlands that extend to areas with multiple systems (Baptist et al., 2019) and cover complete systems (Ysebaert et al., 2003). The generic processes that form and populate estuaries could also be studied by numerical modeling, that increasingly incorporates interactions with eco-engineering species (e.g., Brückner et al., 2020Brückner et al., , 2021 and by landscape experiments that may incorporate eco-engineering plant species (e.g., Kleinhans et al., 2015;Leuven et al., 2018b;Braat et al., 2019;Lokhorst et al., 2019). However, such modeling and experimentation also requires field data for validation. To compare systems and biogeomorphological models objectively and to monitor changes accurately through time, we urgently need detailed data and methods of biogeomorphological classification at the scale of entire estuarine systems. The quantification of system-scale mapping of patterns by data is of immediate practical importance for protection of estuarine habitats worldwide. In the EU this is imposed by law, such as the European Natura2000 legislation, but also for coastalprotection managers from reliable and up-to-date estuarine maps are critical. Generic insights from similar systems would aid habitat protection and restoration, as would insights into what causes the differences in habitat area and heterogeneity between estuaries. Realistic goals for the areal extent and the distribution of estuarine habitats may be based on such trends, and are important to have in view of the competition between societal functions of flood safety, transport and access to ports and ecosystem health. Sustainable management of the estuaries, for instance by dredging and disposal schemes that incorporate ecological goals, requires data of the hydrogeomorphological processes and of habitats. Furthermore, policy documents frequently refer to the importance of gradual water-land transitions on the basis of ecological literature, but there is at present neither quantitative data of such transitions, nor an objective method to collect these.

Ecological Significance of Habitat Areas and Transitions Between Habitats
The highest macrobenthos diversity and biomass is found in intertidal areas (Ysebaert et al., 2003), making them important feeding areas for fish during high tide (Ysebaert et al., 2000).
Tidal marsh creeks are important nursery grounds for some fish and macrobenthic species (Ysebaert et al., 2016). During low tide the intertidal areas rich in macrobenthos are important feeding areas for birds, especially wader birds. Wader birds often forage along the waterlines of tidal flats and require a long feeding period during low water, and longer with lower body weight (Ysebaert et al., 2000). This means that the quantity of low-energy intertidal area is important, as well as the spread of its height distribution relative to the tidal levels. The saltmarsh provides breeding grounds for birds and may serve as resting ground during high water (Ysebaert et al., 2016), showing that animal species move across multiple habitat transitions.
Direct quantitative relations between species abundance and habitat area are simplistic, because species abundance depends on multiple variables such as tidal and seasonal dynamics (Cozzoli et al., 2013). Furthermore, ecological studies also demonstrated the importance of transitions between different habitats (or ecotopes) along the river-estuary transition for species occurrence and abundance. However, quantification of patch variation and gradients for entire landscapes remains a challenge, and many possible metrics with ecological connotations exist (Lausch et al., 2015). It is, regardless of interpretation, clear that species need such transitions in various life stages for feeding and reproduction cycles as discussed below, and that quantification of the transitions is useful.
Relations for individual species are likely to differ between areas, and the use of such species-specific relations for system management would require potentially controversial choices for targeted communities and indicator species by ecologists and policy makers. Instead, we limit ourselves to a minimal review of possible qualitative links between groups of species and the presence of biogeomorphological areas and transitions. The latter are empirically much more accessible and comparable between different systems. This points at the need for objective quantification of transitions between habitats. Habitat diversity and presence of transitions are important for the many species that rely on more than one habitat along estuaries during their life cycle. While most freshwater species do not enter saline waters, many marine species can extend into brackish waters and some fish species need the gradual transition for migration (Ysebaert et al., 2016). The tidal marsh creeks and the transition from the intertidal areas are important as nursery grounds for some fish species (Ysebaert et al., 2016), likely including the pioneer saltmarsh zones as transitions and partially sheltered zones. Possibly, the weakly brackish zone has a smaller number of animal species than the saline or freshwater zones (Whitfield et al., 2012), but the diversity of the primary producers may be larger than elsewhere (Telesh et al., 2013), which is important for macrobenthic and fish species (Herman et al., 1999). Together, this suggests that not only the area of low-energy intertidal flats but also the length of the transitions between intertidal habitats and deeper water are relevant measures for habitat availability.

Objective
The first objective of this paper is to compare trends of biogeomorphological areas along three estuaries and assess whether these are generally comparable or mainly determined by system-specific histories. The second objective is to quantify lengths of transitions between biogeomorphological classes and compare these between the estuaries. To this end, we present a novel image classification method that allows a detailed objective comparison of trends along three estuaries in a similar setting of tidal range and climate, similar shore-connected areas as well as mid-channel bars, but with different forcing histories. Three estuaries in the Netherlands have sufficient data of high quality for automatic image analysis at the required scale. While this low number limits the general applicability to estuaries worldwide, all three have similar boundary conditions, saltmarsh and benthic species and have contrasting historic developments and geologic constraints. For one of the estuaries, the new biomorphological map is compared with a manually created map as a moderate assessment of the accuracy.

Approach
The analysis of habitat indicators along and across estuaries requires consistent data that cover entire systems and generalize the environmental conditions. Not surprisingly, there is very little data published in the scientific literature that cover entire systems. There are only few bathymetries (Leuven et al., 2018c), maps of mud presence (Braat et al., 2017), saltmarsh (Wang and Temmerman, 2013;van Belzen et al., 2017) and macrobenthic species (Cozzoli et al., 2013), mostly in northwestern European estuaries.
Here, we take advantage of our newly developed method of automated object-based aerial image classification for biogeomorphological maps. We applied this to all relevant tidal systems in the Netherlands: the Wadden Sea (in progress) and three estuaries ( Table 1), namely the Eems-Dollard, the Eastern Scheldt and the Western Scheldt (Figure 1, see Table 1 for description of conditions and contrasting histories). For system-scale management purposes and monitoring for the Natura2000 legislation, Rijkswaterstaat has used ecotope maps (Bouma et al., 2005), where "ecotope" is defined as the physical environment for communities, which is akin to biotope and habitat (see Baptist et al., 2019, for discussion). The ecological significance of all ecotopes for the estuaries and the Wadden Sea in the Netherlands is reviewed on the basis of observations and discussed for target species in Bouma et al. (2005), who also acknowledge that boundaries on gradually varying surfaces are artifices on transitional surfaces but are needed for mapping, management and legal purposes (see discussion). Important classes of ecological indicators are saltmarsh and unvegetated tidal flats subdivided in low-dynamic and high-dynamic littoral zones ( Table 2).
Ecotope maps can in part be generated from datasets of bathymetry, salinity, and hydrodynamics. Two methods have been applied to make ecotope maps: one based on numerical modeling of hydrodynamics, and one based on biogeomorphological characterization of the littoral zones from aerial photography. The first method relies mainly on numerical hydrodynamic modeling (Ysebaert et al., 2016;Cozzoli et al., 2017;Baptist et al., 2019) for quantification of conditions in 1 | Key characteristics of the three estuaries: length of channel centerline (also shown in Figure 4, total mapped area (without the water class), spring-tidal range at the mouth, mean and maximum river discharge and general description (see text).  the littoral zone. Such modeling is not only computationally expensive but is also inaccurate for the shallow areas of interest due to the large effects of wind and of chosen hydraulic resistance values, and the spatiotemporal propagation of errors in measured bathymetry and boundary conditions of tides and wind (Duran-Matute et al., 2014). Nevertheless, this is the most straightforward way to determine inundation duration, which is an important variable in habitat suitability and ecotope definition (Bouma et al., 2005;Cozzoli et al., 2013). Determining inundation time in estuaries with progressive tidal waves requires accurate bathymetry, which requires serious surveying time. It also requires many water level measurements and an accurate numerical model that is usually calibrated to measured water levels by adjusting hydraulic resistance spatially. The second method to create ecotope maps is to produce biogeomorphological maps for the intertidal and supratidal areas on the basis of imagery. One difference between biogeomorphological mapping over image interpretation to infer habitat, or modeling of hydrodynamics over bathymetry, is that another kind of information is used about the physical processes that underlie the habitat definitions. Biogeomorphological maps have historically been drawn by hand from aerial photographs based on protocols and checked in the field with a frequency of about once per 2 years for the Western Scheldt, once per 6 years for the Eastern Scheldt, but not at all for the Wadden Sea and Eems-Dollard. The automated mapping method presented in this paper now allows for cost-effective production of biogeomorphological maps in all systems with high-quality imagery. Rijkswaterstaat has a workflow of combining handdrawn biogeomorphological maps with calibrated water level modeling on measured bathymetry to produce ecotope maps for the Western Scheldt every 2 years and Eastern Scheldt every 6 years.
To work with consistent data for the three estuaries, here we used automatically generated biogeomorphological maps, i.e., not ecotope or habitat maps that can be derived from combination of hydrodynamic and biogeomorphological information, that were checked and corrected by experts for reasons of comparability of intertidal and supratidal habitat trends across and along the estuaries (explanation below). The manually produced maps of the Western Scheldt were used in the development stage as verification data for the automated method. We refer to the Supplementary Material for details of the method and a comparison with a manually created map of the Eastern Scheldt as validation.

Study Areas
All three estuaries are situated in the coastal plain of the Netherlands (Figure 1 and Table 1) and formed only late in the Holocene by ingression due to the rising sea level, peat collapse and land subsidence due to dewatering for agriculture (see de Haas et al., 2018, for overview). The climate, the tidal boundary conditions and the river discharge are similar, except that the Eastern Scheldt is strictly no longer an estuary but a tidal basin, since it was disconnected from the Scheldt river. The Scheldt system in particular has been studied extensively and the brief descriptions below are limited to the main characteristics from selected publications. Relevant class abbreviations as also shown in the figures are taken from the original system by Bouma et al. (2005, Appendix 5) whereby ecological significance depends partly on mud content, salinity, pelagic primary production and spatiotemporal variations thereof (e.g., Herman et al., 1999;Ysebaert et al., 2003).
The Eems-Dollard estuary is situated in the Wadden Sea on the Dutch-German border. The lower estuary is not confined or embanked but connected during high water with surrounding tidal basins, and the mouth is confined between two barrier islands. Late-Holocene ingressions led to wide-spread land loss in the area, including a large embayment on the south-western bank of the estuary (Vos and van Kesteren, 2000;de Haas et al., 2018). Gradual embankment and stimulated saltmarsh accretion undid part of the land loss but a large embayment, the Dollard, still exists (Figure 1). The lower estuary has mutually-evasive ebb-or flood-dominated channels separated by tidal bars, where the tidal dominance is enforced by the bends in the system (van Veen, 1950). The mud concentrations in the estuary have increased over the past decades to levels that are considered harmful for primary production (de Jonge, 2000; Winterwerp et al., 2013;van Maren et al., 2015). Causes are not known with certainty but the deepening by several meters of the Ems river upstream of the Dutch-German border is likely to have played a role in the enhanced mud trapping, and both the deepening and the reduced roughness due to increased mud concentrations led to amplification of the tides landward of the Dollard (Winterwerp et al., 2013;van Maren et al., 2015). It is unknown whether this led to differences in ecotopes over time and in comparison to other estuaries in this region, because they were never mapped. The single-thread Ems river in Germany is excluded from our analysis.
The well-known storm-surge barrier of the Dutch Delta Works, built in the 1980's after the great flood of 1953, is located at the mouth of the Eastern Scheldt. Following late-Holocene ingressions, a large area was flooded with some remnant islands remaining. The silting up and damming of minor branches between these islands led to the present shape of the tidal basin with two minor branches to the north (van den Berg et al., 1996;Vos and van Kesteren, 2000;de Haas et al., 2018). These branches were disconnected from the northern basins as part of the Delta Works. Some large tidal bars exist in the estuary (Mulder and Louters, 1994;de Vet et al., 2017). The Eastern Scheldt was connected to the Scheldt river before the Western Scheldt formed in the early Middle Ages (van der Spek, 1997). The present-day lack of fluvial input means that this basin is entirely saline, estuarine circulation is absent and, while it receives only little mud from the seaward boundary, the water is generally clear (compare water color in Figure 1). Damming of the side-branches initially increased the tidal range, but installation of the storm-surge barrier has approximately halved the tidal range at the inlet and reduced the tidal prism and sand input from the North Sea littoral zone much more (Cozzoli et al., 2013). The fact that the tides and dynamic salinity gradient were not simply entirely removed by a dam, as in other systems in the Netherlands, was the result of a compromise between the aquaculture sector and the nature conservation organizations.
The Western Scheldt is situated in the south of the Dutch North Sea close to the Belgian border. As in the Eastern Scheldt case, the landscape in the Middle Ages was much more fragmented with large branches that originally provided nautical access to the late medieval ports of Middelburg, Gent, and Brugge but gradually silted up and were dammed (van den Berg et al., 1996;van der Spek, 1997). In 1570 a large embanked polder, Saeftinghe, was flooded as a defense in a war, and since then silted up to supratidal levels to form a large saltmarsh that was only partly reembanked (Wang and Temmerman, 2013;Jongepier et al., 2015). The mainport of Antwerp necessitates direct access to the North Sea, and a storm surge barrier was rejected after a very long societal and political debate. Tidal channels were deepened for ships that gradually increased in draft over the past decennia. A number of large tidal bars formed over time with a main, ebb-dominated channel and a number of flooddominated channels (Wang et al., 2002;Hibma et al., 2003;Hiatt et al., 2020); a pattern enhanced by large-scale channel deepening (van Dijk et al., 2019). Meanwhile the tidal amplitude increased landward due to the channel fairway deepening (van der Spek, 1997; Wang et al., 2002). Small-scale channels cross-cut the bars as part of the general dynamics (van Veen, 1950;Swinkels et al., 2009), but historic analyses suggest that this reduced over the past half century (Cleveringa and Taal, 2015). The Scheldt river input causes a salinity gradient with the mesohaline zone around the Dutch-Belgian border and the resulting estuarine turbidity maximum close to Antwerp (Winterwerp et al., 2013). In order to mitigate the increasing flood risk, controlled inundation is allowed in designated areas, again after considerable societal and international political debate. The single-thread Scheldt river in Flanders is excluded from our analysis.
The tidal bars with their intertidal flats showed opposite developments in the Western Scheldt and Eastern Scheldt over the past decades (de Vet et al., 2017). The bars in the Western Scheldt generally increased in height since fairway deepening was commenced, and as the system was confined in width, this caused general steepening of the bar margins (van Dijk et al., 2019). In contrast, the bars in the Eastern Scheldt lowered since the storm surge barrier was installed, mainly due to the relative increase in wave attack due to the tidal range reduction in combination with the excess accommodation space in the tidal channels, of which the dimensions were related to the much larger tidal prism preceding the Delta Works. As a result, the bar margins became less steep (de Vet et al., 2017).
All three systems are Natura2000 sites designated under the EU Birds Directive and the EU Habitats Directive, meaning that tidal flats and saltmarsh are important habitats to be monitored (e.g., Meire et al., 2005). The ecological consequences of human interference have largely been negative. The Delta Works have started a gradual decline of the surface area and height of the tidal flats and have limited rejuvenation the salt marshes in the Eastern Scheldt, with a low but steady reduction of bird habitat as a consequence. For the Western Scheldt and Eems-Dollard systems changes in tidal dynamics, salinity and mud concentration as a result of the ongoing morphological response to dredging of fairways in combination with past reductions of intertidal area by land reclamations have resulted in the loss of valuable habitats and a reduction of the quality of these habitats (de Jonge, 2000;Meire et al., 2005;Ysebaert et al., 2016).

Aerial Photographs and Image Classification
The biogeomorphological maps were based on high-quality digital aerial photographs with near-infrared, red and green bands and a ground resolution of 0.25 m. This resolution is much higher than obtained in satellite imagery, and is necessary because the morphological phenomena of interest on which the classification is based are then resolved, or nearly resolved. While Unmanned Airborne Vehicles would obtain an even higher resolution, the spatial extent possible with the aerial photographs is much larger, which is necessary given the focus on entire systems and fluvial-tidal river gradients. We found only the most recent sets to be of sufficient quality for digital classification, in part because of technological progress and in part because only these met the strict requirements of having been collected with a low tidal water level range, in the same season, and with uniform lighting and weather conditions to avoid large differences in reflection intensity and color, and avoiding direct reflection of the Sun.
After extensive testing an Object-Based Image Analysis ruleset was developed, of which a brief summary is provided here (see Supplementary Material for more information). Tiles of merged false-color aerial photographs were segmented multiple times, at different levels of scale, to yield a few tens of thousands of objects with a minimum degree of spectral heterogeneity. The smaller-scale objects were nested in the larger-scale objects. We developed the ruleset in eCognition to classify the objects into seven classes ( Table 2). The classification was based on object shape and context (neighboring objects), relations between the nested objects, on average brightness and on a ratio between near-infrared and red reflection: the Normalized Difference Vegetation Index (NDVI). Distinguishing the mud content of the surface proved challenging for the entire system at this scale, despite earlier success with remote sensing imagery (van der Wal et al., 2010), and was not conducted. Hard stone surfaces, which in these systems all have an anthropogenic origin, were copied from digitized maps of Rijkwaterstaat of the Ministry of Public Works.
The ruleset is summarized as follows (see Supplementary Material for more information). Water was classified by a maximum NDVI. Saltmarsh was classified by a minimum NDVI. Saltmarsh objects smaller than 400 m 2 were dilated and eroded, which removes vegetation patches too small for the mapping scale. Merged saltmarsh objects larger than 400 m 2 were then classified as pioneer vegetation, which is less densely vegetated than open salt marsh. Saltmarsh creek objects were recognized from their unvegetated appearance, elongated shapes and vegetated neighbor objects. The remaining objects were unvegetated tidal flats, for which the two levels of segmentation were used. Mega-ripples were distinguished by variation of brightness and NDVI of sub-objects nested in the larger objects. The final discrimination between high-energetic and low-energetic tidal flat was made by combinations of non-overlapping thresholds of brightness and NDVI, where undecided objects between the thresholds were classified from the surrounding context.
The automated classification based on the ruleset produced a similar map as the visual method (example of the automated map of Eems-Dollard in Figure 2). The spatial detail at which surface areas were classified is illustrated for saltmarsh and tidal flats in Figure 3.
The automation reduced the map-production time considerably. For present conditions and an entire estuary, the total computation time is 24 h. In contrast, manual mapping cost several weeks, particularly because of the need to draw all map boundaries. An expert check of the digital map was fast because the boundaries were found to be accurately captured by the objects, even if objects were not correctly classified. The accuracy of the boundaries is partly determined by the scale number in the object-based image classification: earlier attempts with higher scale numbers occasionally resulted in clearly mislocated boundaries.
The automated classification and the expert-corrected maps of the Wadden Sea systems and the Eastern Scheldt have 90% or more correspondence and about 85% in the Western Scheldt. For saltmarsh areas, including partially covered areas, the agreement is about 95%, where differences arise due to objectified thresholds for partial cover and the detail in the boundaries in the automated  Table 2 for explanation. Full image resolution maps for all three systems are provided in Supplementary Material. method. Most of the disagreement between the maps occurs in the bare surface classes. We tested whether the use of bed slope and bed-slope variation calculated from bathymetry improved the classification of the bare surfaces. Although this increased the correspondence of mega-ripple area on both maps, it added low-energy and high-energy plane shoal areas to the mega-ripple areas that were in disagreement with the manually drawn map. Visual inspection shows that the boundaries between tidal flat and mega-ripples are sometimes poorly recognized in the automated method. The expert-corrected spatial data are made available online and in this paper we proceeded with the corrected data.

Map Analysis
The (corrected) classified objects were used in the following analyses, the information for which is included in the online data. To study spatial trends along estuaries we drew narrow segments across the estuaries and calculated the absolute fraction of the surface area covered by each legend class from the classified objects. These were plotted as absolute widths (surface area divided by the width of the segments) and as fractions of the surface area to visualize their absolute and relative importance along the estuary.
All three estuaries are characterized by the presence of large mid-channel bars. To test the hypothesis that such bars are more dynamic environments than the shore-connected tidal flats on the estuary margins that are partly protected from tides and waves, we also plot trends along the estuaries for the mid-channel intertidal areas and the shore-connected areas.
A measure for spatial habitat variability is the length of transitions between biogeomorphological classes. After merging objects of the same classes, the lengths of the borders between areas of different classes were calculated for the entire estuaries and for the mid-channel intertidal areas and the shore-connected areas. The resulting absolute lengths obviously depend on the scale of the Object-Based Image Analysis (OBIA), on the resolution of the images and on the actual spatial variation. The two method-related factors are the same for all three estuaries in this study. The pixel size of the images (0.25 m) used here is at a similar scale as the birds and fish, so that the analysis of spatial heterogeneities and lengths of transitions is done on a similar scale as that of the locomotion paths of birds and fish. This suggests that the quantitative measures presented here are likely relevant for the organisms.
To contextualize the transition lengths, we normalized the more abundant, measured transition lengths by a length that is characteristic for the total length of the mapped area along each estuary. This characteristic length was calculated from the total mapped surface area of the estuary (without the water) divided by the mouth width (Table 1). We interpret the normalized transition length as the minimum distance that a transition between saltmarsh and tidal flats, and between tidal flats and water, would have in a highly idealized and regular estuary. The more irregular the mapped transitions between habitats are, the larger the normalized transition length. While the irregularity and length of transitions between tidal flat classes and between saltmarsh classes is at least partly determined by the definitions of the classes and the OBIA technique, the transition lengths between saltmarsh and tidal flat classes are quite likely robust due to the relatively straightforward recognition of vegetation on imagery with the near-infrared band. Furthermore, the image resolution and OBIA ruleset were identical in the three study areas.

Trends Along Estuaries
The relative area covered by the biogeomorphological classes show two trends that are similar between the three estuaries. First, there is much variation along the estuaries in the intertidal and supratidal width fraction (Figure 4) mainly related to the presence of large bar complexes (Figure 1). In all cases, the subtidal takes up the largest area, the tidal flats are second and the vegetated areas are smaller (Figure 5).
The second similarity is that the proportions of tidal flat and vegetated area reduce in seaward direction (Figure 4) but the absolute areas remain about the same (Figure 5). In the absolute sense, the most upstream reaches of the Western Scheldt and Eems-Dollard are in places wider than the estuary mouth, and in all three estuaries the landward ends have a disproportionally large area above the subtidal level.
There are three striking dissimilarities between the three estuaries. First, there is proportionally much more high-energetic tidal flat and megaripple fields in the Western Scheldt than in the other two estuaries (Figure 4). Second, the Eastern Scheldt and Western Scheldt contain a small but important fraction of hard substrates of natural origin, mainly compacted, early-Holocene peat layers (Figure 4).
The third dissimilarity is that much of the upstream area is vegetated in the Western Scheldt (Saeftinghe salt marsh), while in the other two estuaries low-energy tidal flat dominates the upstream area (Figure 4). In terms of absolute width, the vegetated area in the upstream reach of the Western Scheldt is almost as wide as the subtidal channel and the unvegetated flats together, as opposed to everywhere else in the estuary (Figure 5). Generally, the Western Scheldt has more vegetation than the other two estuaries. The Eastern Scheldt has two side branches that both have more vegetation than the main branch (Figure 6).
Two specific features in the Eems-Dollard are related to its position in the Wadden Sea and the somewhat arbitrary choice what to incorporate in the Eems-Dollard study area, because its boundaries within the Wadden Sea are hard to define and we used the administrative boundary of the area for simplicity. The large vegetated patch near the mouth of the Eems-Dollard is located on a barrier island in the Wadden Sea that is arguably not part of the estuary. The large widening of subtidal width in the Eems-Dollard between 5-10 km from the mouth is the historic connection to the eastern branch of the Eems estuary ( Figure 5).
Unexpectedly, there is not proportionally more vegetation in the more sheltered shore-connected areas than on the tidal midchannel bars, despite the fact that the latter are surrounded by channels and are less sheltered (Figure 7). The most obvious trend is a general reduction of shore-connected area and an increase of mid-channel bars associated to the seaward widening of all estuaries. However, there is more high-energetic and megaripple area on the mid-channel bars of all three estuaries in general and particularly in the Western Scheldt, except in the most upstream, single-thread reaches. Most of the hard substrates in the Eastern Scheldt are located on the shore-connected areas.

Transitions Between Biogeomorphological Classes
The sum of lengths of all the transitions is well above 2,000 km (Figure 8), which is, not surprisingly for a patchy landscape, about an order of magnitude larger than the circumference of each entire estuary. The analysis shows that a minority of all possible transitions dominates in all three systems: from saltmarsh (S1a) to pioneer vegetation (S2) and saltmarsh creek FIGURE 4 | Distribution of biogeomorphological classes as relative width along the Western Scheldt (top), Eems-Dollard (middle), and Eastern Scheldt (bottom). Cumulative estuary width, the sum of all classes including water, is plotted as the dashed black line for the right vertical axis. Distance is measured along a smooth centerline between the banks from the narrowest cross-section in the seaward mouth. See Table 2 for explanation of the legend. The data underlying these graphs are provided in Supplementary Material.
(S3a), from pioneer vegetation (S2) to low-energy tidal flat (P1) and from low-energy tidal flat (P1) to high-energy tidal flat (P2c), anthropogenic hard substrate (H2) and water. As expected, the saltmarsh classes occur frequently adjacently and saltmarsh also borders frequently on low-energy tidal flats. Furthermore, low-energy tidal flats border frequently on high-energy tidal flats, anthropogenic hard substrates such as bank protections and water.
On the other hand, some general differences between the estuaries emerge (also see Figure 9). The Western Scheldt has much more transition length between mega-ripples (P2b), highenergy tidal flat (P2c) and water than the other two, because the mega-ripples are more abundant in this system. These are also the classes for which the discrimination from low-energy tidal flats was more difficult than in the other two estuaries (see section 2). On the other hand, the Eastern Scheldt has much less transition length from saltmarsh (S1a) to low-energy tidal flat (P1) and from high-energy tidal flat (P2c) to water. In general, the Eems-Dollard can be said to take up an intermediate position between the Western Scheldt and Eastern Scheldt in absolute transition FIGURE 5 | Distribution of biogeomorphological classes divided by vegetated and unvegetated parts, plotted as absolute width along the Western Scheldt (top), Eems-Dollard (middle), and Eastern Scheldt (bottom). Positive widths refer to subtidal bed surface and to all vegetated classes, negative widths refer to unvegetated intertidal flat classes and the sum of absolute widths is the cumulative estuary width excluding water. Distance is measured along a centerline from the narrowest cross-section in the seaward mouth. For explanation of legend see Table 2. lengths between saltmarsh (S1a) and high-energy tidal flat (P2c), and between low-energy tidal flat (P1) and mega-ripples (P2b).
The lengths of transitions are different for shore-connected areas and mid-channel bars and shoals (Figure 9), which is obvious from the fact that the matrices are not perfectly mirrored on the diagonal. Main similarities between the three estuaries are that tidal flats (P1, P2b, P2c) border more frequently on water for the more exposed mid-channel bars than for the shore-connected areas, which supports the hypothesis that shore-connected areas are less dynamic than mid-channel areas. Furthermore, vegetated surfaces border on unvegetated surfaces about equally frequently on the bars and shore-connected areas. The more frequent occurrence on the shore-connected areas in all estuaries of lowenergy tidal flats bordering on anthropogenic hard substrates FIGURE 6 | Distribution of biogeomorphological classes separately plotted as cumulative width for the main branch (top) and two side branches of the Eastern Scheldt, one with a SW-NE orientation (middle), and the small W-E channel (bottom, middle channel in Figure 1). The sum of the branches is plotted in Figure 5 (middle). Distance is measured along a centerline from the narrowest cross-section in the seaward mouth. For explanation of legend see Table 2.
(H2), which includes bank protection works and infrastructure, is unsurprising.
The most abundant transition lengths were normalized by a characteristic length of each estuary for comparison between the three estuaries ( Table 3). In all estuaries the transitions between saltmarsh classes and low-energy tidal flat are an order of magnitude longer than those with the high-energy tidal flats. The relative transition lengths between low-energy tidal flat and saltmarsh classes differ between the estuaries: the Western Scheldt has relatively more transition length to saltmarsh than to pioneer vegetation, while the Eastern Scheldt has relatively much more transition length to pioneer vegetation. The Eems-Dollard has the lowest relative transition lengths of all three estuaries.
There are striking differences between the estuaries that are consistent with the differences observed before. The Eastern Scheldt has no significant vegetated area on the bars, in contrast to the other two estuaries, and therefore has no FIGURE 7 | Distribution of biogeomorphological classes separated for shore-connected and mid-channel areas as absolute width along the Western Scheldt (top), Eems-Dollard (middle), and Eastern Scheldt (bottom). Positive widths refer to shore-connected tidal flats and saltmarsh and negative widths refer to the bars and shoals. Distance is measured along a centerline from the narrowest cross-section in the seaward mouth. For explanation of legend see Table 2. significant transition length between vegetated and unvegetated habitats. The Western Scheldt has more shore-connected megaripples (P2b) and high-energy tidal flat (P2c) bordering on water, while the Eastern Scheldt has more low-energy tidal flat (P1) bordering on water and the Eems-Dollard takes an intermediate position. On the other hand, the Eastern Scheldt and Eems-Dollard have much more transitions from low-energy tidal flat (P1) to water on the bars than the Western Scheldt.

DISCUSSION
The main outcome is that the three estuaries have fairly similar biogeomorphological characteristics with predominantly subtidal areas, here unmapped, and a larger littoral zone of unvegetated tidal flats, than pioneer saltmarsh and supralittoral saltmarsh. The tidal flat width varies along the estuaries, which was already explained by formation of mid-channel bars and shore-connected embayments in zones wider than  Table 2. The color is linearly scaled with the numbers.
Frontiers in Marine Science | www.frontiersin.org 13 October 2021 | Volume 8 | Article 742448 FIGURE 9 | Fraction of transition length between biogeomorphological surface classes for shore-connected areas (red) and islands and shoals (blue) (Continued) FIGURE 9 | for the Western Scheldt (top), Eems-Dollard (middle), and Eastern Scheldt (bottom). For explanation of classes see Table 2. The colors are linearly scaled with the numbers. The sum of each entire matrix is 100%.
a monotonically converging minimum estuary width (Leuven et al., 2018b,c). The quality of the classification is as good as expected from similar applications (e.g., Owers et al., 2016).

Interpretation in the Context of Contrasting Histories
The spatial variation of saltmarsh area is in all three cases most likely related to the processes of estuary formation by ingression on an inherited, embanked landscape, and by the land reclamation history. Historic and geologic reconstructions showed that parts of the presently embanked area were saltmarsh, and parts are former side-branches that silted up or were closed off (de Haas et al., 2018). The precise situations differ between estuaries: the Western and Eastern Scheldt had extensive side branches (Mulder and Louters, 1994;van der Spek, 1997). The majority of the upstream saltmarsh of the Western Scheldt developed over the past three centuries in Saeftinghe, a former polder flooded as a defense measure during a war. And the upstream saltmarsh in the Eems-Dollard formed in the Dollard, an inland bay that was a remnant of the ingression. Much of the original bay was reclaimed after saltmarsh formation intensified by human-built wave breakers with the intention to gain land, which was mainly successful in the predominantly downwind direction (Vos and van Kesteren, 2000;de Haas et al., 2018). It is not entirely clear, however, why Saeftinghe in the Western Scheldt silted up to high supratidal levels without human-built wave breakers, while the Dollard ceased to expand and fill. Spartina anglica, a highly effective species for humanaided salt marsh formation, was introduced and spread in the 1920s. There are also local dynamics in the saltmarsh on shorter timescales, such as quasiperiodic retreat and expansion (Wang and Temmerman, 2013;van Belzen et al., 2017), but these do not change the general trends. The patch of saltmarsh near the mouth of the Eems-Dollard is more associated to the barrier island than to the estuary. This leads to the interpretation that the large-scale pattern of landward increasing saltmarsh extent in all three estuaries is mainly caused by the initial landscape conditions prior to the ingressions and by the land reclamation history of the past centuries. Human interference on the much shorter timescale of the past half century was also expected to have affected the littoral biogeomorphological patterns in all three estuaries. In the Eastern Scheldt, the installation of the storm surge barrier in the mouth instantaneously reduced the tides, which not only reduced the dynamics in the sea arm but also led to channel sedimentation, bar erosion and reduction of the littoral zone (Mulder and Louters, 1994). Indeed, dredged sand has now been dumped on the largest, western bar to attempt stalling of the loss of littoral habitat due to the reduced tidal energy in the adjacent channels. However, compared to the other two systems, the largescale pattern is not yet very different. The Eastern Scheldt has Measured transition lengths were divided by a characteristic length calculated from the total mapped surface area of the estuary (without the water) divided by the mouth width (Table 1). See Table 2 for explanation of ecotope labels.
a much larger area of hard substrates than the other estuaries, which is likely partly due to the antecedent landscape conditions and partly due to loss of the sediment cover from these hard substrates to the subtidal zones. The intertidal areas of the Eems-Dollard are not notably different from those in the other systems, but note that our method cannot yet distinguish between sandy and muddy areas. Particularly the Dollard embayment has a very high mud content (van de Lageweg et al., 2018). The two most striking differences between the systems are the presence of saltmarsh on the mid-channel tidal bars and the large high-energy area in the Western Scheldt. The formation of saltmarsh on the bars is likely an effect of dredging, fairway deepening and disposal of sediment on bar margins. The deepening enhanced the tidal amplitude throughout the estuary (van der Spek, 1997;Wang et al., 2002Wang et al., , 2015. While adaptation of the entire streamwise profile to the changed tidal properties and continuously disturbed sediment balance may take centuries (van der Wegen and Roelvink, 2012), the adaptation of the local "transverse" profiles and morphology were faster. The enhanced tidal currents scoured the channels further and transported sediment onto the bars, in addition to the disposal of dredged material on bar margins (van Dijk et al., 2021). Around the turn of the century the mid-channel bars became so high that the reduced inundation times allowed saltmarsh to settle and expand. Unlike the shore-connected intertidal flats and saltmarsh, this portion of the saltmarsh area is unrelated to the inherited landscape and land reclamation histories but is directly related to altered sediment dynamics. The reasons for the much larger high-energy area in the Western Scheldt are poorly understood, but the high proportion of mega-ripple fields predates the human interference (Cleveringa and Taal, 2015).

Significance for Potential Habitat Areas and Transitions
The review and the mapping exercise in Bouma et al. (2005) and in this paper emphasize the quantification of areas of different biogeomorphological units as well as their transitions. While no specific habitat is strongly preferable from an ecological perspective, a mosaic of habitats and transitions may be preferred over monotonous landscapes by animal species from several taxonomic groups. Ecological mapping has been argued to fall into either the category of patch delineation or the category of landscape gradient mapping (Lausch et al., 2015). The present mapping approach clearly falls into the category of quantification of habitat area by patch delineation, but the new measure for length of transition between habitat classes also honors the ecologically important idea of landscape gradients.
The low-energy intertidal area appears more important than the megarippled and high-energy tidal flat areas for abundance and biodiversity of macrobenthic species and their predators (Ysebaert et al., 2003(Ysebaert et al., , 2016. The sediment dynamics of the high-energy tidal flats restrict the number of benthic species that can successfully survive in and on these flats, and therefore also limits the foraging time for waders during the emergence period. Furthermore, it suggests that the large proportion of highenergy vs. low-energy intertidal area in the Western Scheldt is less desirable from the perspective of preservation in view of reducing habitat areas regionally and globally. Regardless of the abundance of high-energy habitats in the Western Scheldt, the ecologically important habitats of saltmarshes, with relatively large lengths of transitions between low-energy tidal flats and saltmarsh creeks, occur more frequently than along the Eastern Scheldt and Eems-Dollard. The Eastern Scheldt has limited saltmarsh and reduced intertidal flat area due to the effects of the storm surge barrier (de Vet et al., 2017). The Eems-Dollard has limited saltmarsh with creeks because of human interference during saltmarsh accretion for embankment purposes. Furthermore, much of the lower estuary here is part of the Wadden Sea and therefore unconfined by shores where saltmarsh could have formed.
The presence of the large mid-channel bars, as opposed to single-channel estuaries, increases the intertidal area, the length of water-tidal flat transitions and the transitions from low-energy tidal flats to saltmarsh creeks in case saltmarsh develops on the bars. While the shore-connected tidal flats and saltmarshes along most of the systems may be more sheltered than their midchannel counterparts, the total length of transitions on the sparse bars is as important numerically and, for the saltmarsh creeks, more important. More generally, this shows that multi-channel estuaries with bars have potentially much richer and more extensive habitats than estuaries that were confined to a single channel system for shipping and land-reclamation purposes. Additionally, access to mid-channel bars by land predators, including domesticated animals, is likely much lower than to the shore-connected areas.
The lengths of the most abundant mapped transitions between tidal flats and saltmarshes were compared after normalization (see Table 3). As expected, the normalized transition lengths are high for the transition between low-energy tidal flats and natural pioneer vegetation and saltmarsh creeks, but much less so in the Eems-Dollard than in the other two systems. This is consistent with the higher human influence on saltmarsh formation in the Eems-Dollard, which renders the transitions between bare and vegetated surfaces much more regular and relatively shorter. On the other hand, the tidal flat transition to fully developed saltmarsh is especially long in the Western Scheldt, which fits with the large saltmarsh area in that system. The transition from high-energy tidal flats to pioneer vegetation is similar in the three estuaries after the normalization, despite the much higher abundance of high-energy tidal flats in the Western Scheldt. Overall, the relative length of transitions between low-energy tidal flat and the collective saltmarsh classes is about 4-5 times larger in the Western Scheldt and Eastern Scheldt than in the Eems-Dollard. In view of the policy that a range of habitats and spatial heterogeneity is required for ecologically meaningful transitions between the ecotope classes (Wolter et al., 2016;Ysebaert et al., 2016), this could suggest that the Eems-Dollard has fewer habitat transitions than possible, as observed in the other two estuaries.

Perspectives
The automated mapping of biogeomorphology from aerial photography is in principle repeatable for multiple years and for other areas, provided that aerial photographs of sufficient quality are available. The photos are ideally acquired on a single day or limited number of days, during low tide and with thin clouds higher in the atmosphere. We ordered orthomosaics with high geometric accuracy and processed them per acquisition day. The model starts with a calibration on Brightness and NDVI, which proved to be the main factors responsible for differences between the photos. The subsequent rules needed minor calibration. For the distinction between mud and sand, brightness appeared to have insufficient information: while mud often has a different moisture content and hence brightness than sand, more factors affect the brightness threshold than we have control over for the system-wide mapping, such as time since emergence and subtle lighting and bedform differences. Future extension of the method to distinguish between sediment types will likely require new techniques. Further extension to distinguish between vegetation types is also an interesting improvement to pursue (Owers et al., 2016). The present mapping will need to be complemented by further analyses and data to obtain habitat maps, or ecotope maps, that include hydrodynamics and salinity dynamics (e.g., Baptist et al., 2019). The natural biogeomorphodynamics of the systems are poorly quantified but significant within the seasons (e.g., Swinkels et al., 2009;van Dijk et al., 2019;Brückner et al., 2021), at least an annual repetition of the mapping is highly recommended. A comparison between sample data of surface and sediment properties and species with the biogeomorphological map could facilitate interpolation for mapping of the substrate and species data.
The high economic pressures in the three studied systems led to an accelerated deterioration of the ecological quality over the past decades, while the ecological importance for the rivers, the North Sea and the Eastern Atlantic Flyway of migratory birds is immeasurable. It is convenient for purposes of management and legislation to have a mapping method with objective and quantitative characteristics, because, alternatively, uncertainty about such numbers could lead to conflicts and legal disputes about conflicting economic and ecological interests (van Buuren et al., 2010). While the ecological significance of the biogeomorphological units, and the transitions between them (Soininen et al., 2015), need further quantification by comparison to ecological datasets and numerical, calibrated flow modeling (e.g., de Jonge, 2000;Ysebaert et al., 2003Ysebaert et al., , 2016Cozzoli et al., 2017;Brückner et al., 2020Brückner et al., , 2021, qualitatively it is certain that they are important. The observed degradation of ecosystems leaves no doubt that effects of human interference urgently needs to be reduced and its effects reversed in order to stop the rapid decline of biodiversity (also see Crooks and Turner, 1999;Ysebaert et al., 2016). Habitat protection and expansion by "depoldering, " or managed realignment, are essential to maintain and restore healthy ecosystems (e.g., Ysebaert et al., 2016). Quantification of areas of habitats, or ecotopes, has already been an efficient monitoring instrument in the Western Scheldt in this respect and the automated method employed here reduces the production time and costs considerably and potentially provides an important spatial dataset for analyses and, possibly, upscaling of ecological datasets.

CONCLUSIONS
A comparison of biogeomorphological maps of three estuaries with similar boundary conditions showed that saltmarsh decreases weakly in width proportion along the estuaries while tidal flats have similar area along the estuary. The subtidal areas have the largest surface area in all three systems, followed by tidal flats and salt marsh. The variations in their relative contribution along the estuaries are due to presence of mid-channel bars and shore-connected embayments.
However, the differences are more striking and are mainly due to the different histories of estuary formation and human interference. The estuary outlines are determined by their ingression and the subsequent land reclamations and embankments over the past centuries. Consequently, the heads of all three estuaries have disproportionately more area above subtidal, where one estuary has much more saltmarsh and the other has more low-dynamic tidal flat. Likewise, side-branches and embayments are remnants of past larger branches and their planforms are the result of the embankment history. The minor importance of saltmarsh flanking the estuaries is due to the land reclamation that turned most of the surrounding marshes into agricultural areas.
Human interference over the past decades strongly affected the areas of, and transitions between, biogeomorphological units as is illustrated by the contrasts between the three systems. The Western Scheldt has longer transitions between low-energy tidal flats and the saltmarsh creeks, that are nursery grounds, not only because of the saltmarsh developed on an abandoned polder, but also because the increasing tidal amplitude along the estuary that led to mid-channel bar accretion and tidal marsh settling. The Eastern Scheldt, on the other hand, had reduced tidal energy that leads to reduction of low-energy tidal flats. The Eems-Dollard, finally, has much less saltmarsh-creek area and proportionally smaller lengths of low-energy intertidal area to saltmarsh creek transitions due to the human-enhanced saltmarsh accretion for land reclamation purposes.