On the Impact of Salt Marsh Pioneer Species‐Assemblages on the Emergence of Intertidal Channel Networks

Tidal marshes play an important role in climate change mitigation through natural coastal protection. The effectiveness of the natural coastal defense by tidal marshes is closely related to their channel network which is in turn greatly influenced by their vegetation cover and shape. Previous research suggests a dual effect of vegetation on marsh topography; stabilizing sediment on the one hand versus promoting erosion and channel incision on the other hand. This study links these effects to different vegetation species, Salicornia procumbens, Spartina anglica, and Puccinellia maritima (further referred to as Salicornia, Spartina, and Puccinellia), by means of a coupled bio‐hydromorphodynamic modeling study. Single species, species‐assemblages, and species shifts were studied, incorporating both species‐specific physical plant properties and spatiotemporal growth strategies. The results indicate the influence of vegetation on the marsh topography to be highly species‐dependent, but also of a very complex nature. Both the presence of Spartina and Puccinellia resulted in significant channel development, whereas Salicornia did not induce topographic change. The combination of several species promoted or reduced channel development depending on the included species. Species‐shifts linked with climatic changes resulted in increased erosion of the existing channel network potentially reducing the protective capacity of the marsh.


Introduction
Coastal ecosystems are some of the most valued and exploited natural systems around the world (Costanza et al., 1997). Among them, tidal salt marshes are of specific relevance for their role in natural coastal protection . For instance, salt marshes can attenuate storm surges and waves through friction induced by their aboveground biomass. Besides, they provide extra storage area for water and thereby protect the hinterland from flooding in converging estuaries (Leonardi et al., 2018;Temmerman et al., 2013). In respect to global change, tidal salt marshes are capable of following sea level rise due to enhanced sediment deposition during the extended flooding periods; thus, building natural levees (Kirwan & Megonigal, 2013;Morris et al., 2002). Accordingly, the coastal protection function of tidal marshes is not only crucial for present coastal populations, but also for generations to come (Leonardi et al., 2018).
Tidal marshes consist of vegetated marsh platforms intersected by tidal channels (Allen, 2000). Tidal channels, being the main conduit for the exchange of water, sediment, and organisms to the adjacent open water body, provide crucial habitat structure and heterogeneity for tidal marshes. Channels determine aquatic and terrestrial community composition (e.g., Visintainer et al., 2006), hydrodynamics (e.g., Temmerman et al., 2005), and sediment budgets (e.g., French & Spencer, 1993). Moreover, the effectiveness of the coastal defense function exerted by tidal marshes is closely related to their channel network (Leonardi et al., 2018). For instance, channel dimensions directly relate to landward flood propagation (e.g., Leonardi et al., 2018;Stark et al., 2015), while channel density and extent influence sediment distribution and the ability of the marsh to keep up with sea level rise (Kirwan & Megonigal, 2013;Leonardi et al., 2018;Sanderson et al., 2000). Although the channel network affects the distribution of vegetation in the marsh, vegetation colonization also influences channel formation. Locally, roots have a stabilizing effect on sediment and aboveground biomass increases flow resistance, resulting in increased deposition and stabilization of already existing channel banks, further referred to as short-scale positive feedback (Coco et al., 2013;D'Alpaos et al., 2006). On larger scales (vegetation patch-scale), vegetation has the opposite effect on flow as flow deviation around dense vegetation patches leads to flow acceleration, which results in increased erosion between patches that can lead to channel incision, referred to as large-scale negative feedback (Schwarz et al., 2014;Temmerman et al., 2007). The balance between these mechanisms, also referred to as scaledependent feedbacks, is closely related to ecosystem properties such as sediment properties, hydrodynamics, and plant characteristics (Schwarz et al., 2018;Van Wesenbeeck et al., 2008).
Previously, interactions between vegetation, flow patterns, sediment dynamics, and thus tidal channel formation have been attributed to physical plant properties mainly representing a generic Spartina-type salt marsh species. For instance, Temmerman et al. (2007) showed how temporal increase in plant density influenced scale-dependent feedback strength and consequently increased channel density using Spartina anglica as a blueprint. Other studies modeled salt marsh morphology based on Spartina alterniflora as a model species (Best et al., 2018;D'Alpaos et al., 2007;Kirwan & Megonigal, 2013;Morris et al., 2002;Rodríguez et al., 2017). However, physical plant properties and spatiotemporal growth strategies vary widely among different salt marsh pioneer species, rendering the general Spartina-type case an oversimplification. Species-specific physical properties, such as stem height, stem density, or stem diameter, result in differing magnitudes of patch adjacent velocity acceleration . Species-specific spatiotemporal growth strategies result in differences in vegetation cover relative to hydrodynamics stresses over time (Schwarz et al., 2018). Accordingly, tidal marshes covering a stress gradient from low-intertidal to supratidal are rarely occupied by a single species (e.g., Allen, 2000). For instance, the common salt marsh grass Puccinellia maritima is relatively sensitive to long periods of inundation and therefore restricted to higher marsh elevations while the halophyte Salicornia procumbens is particularly abundant in the lower intertidal zone. Previous studies have shown that this balance between stress-tolerance and competition across the intertidal stress gradient influences sedimentation and erosion patterns (Bockelmann et al., 2002;Erchinger, 1985;Pennings & Callaway, 1992). This leaves the question how this relative complexity of the real situation influences the emergence and shape of intertidal channels in the presence of different species and multispecies assemblages. In addition, due to global climate change topographical shifts in plant distribution are expected (e.g., Kelly & Goulden, 2008;Lenoir et al., 2008). For example, Spartina is likely to benefit from rising temperatures as this enables earlier germination and results in a head start compared to other species (Gray et al., 1991;Gray & Mogg, 2001;Loebl et al., 2006). On the other hand, increased drought will benefit the salt-tolerant succulent Salicornia at the expense of other common salt marsh grasses (Strain et al., 2017). Such topographical species shifts will have a substantial effect on vegetation composition, channel development, and provided ecosystem functions of tidal marshes (Gray & Mogg, 2001;Strain et al., 2017).
This study investigates the impact of three dominant marsh species of NW Europe (and combinations thereof) on tidal channel formation: Salicornia procumbens, Spartina anglica, and Puccinellia maritima . The objective of this study is threefold; (1) investigate the effects of different speciesspecific vegetation-traits on the initial development of channel drainage patterns, considering speciesspecific physical plant properties as well as spatial and temporal variation in growth strategies (i.e., life history strategies); (2) elucidate the impact of species assemblages (different combinations of salt marsh species) on channel development; (3) investigate the impact of shifts in dominant salt marsh species on channel and vegetation development.

Materials and Methods
The effect of varying species traits on channel development is investigated by means of coupling between a vegetation model operated in MATLAB and the hydro-morphodynamic model (Delft3D) further referred to as bio-morphodynamic model (Hydraulics, 2006

Hydro-Morphodynamic Model
Delft3D is a widely used open-source model environment able to calculate flow, sediment transport, and morphological change and has been applied and validated in various coastal and estuarine environments (e.g., Braat et al., 2018;Fagherazzi et al., 2014;Schwarz et al., 2014). Delft3D simulates flow by solving the conservation of mass and conservation of momentum equations. These equations describe the temporal velocity variations in the x-y plane as a function of advection, friction, eddy diffusivity, water depth, and streamline curvature (Lesser et al., 2004). The FLOW-module of Delft3D was used to simulate unsteady flow and sediment transport on a rectangular grid in 2DH (depth-averaged) driven by a harmonic tidal forcing at the open boundary (Hydraulics, 2006).
For computational efficiency and to ensure tidal propagation independent of lateral boundaries, modeling is performed using two grid resolutions, a small fine grid (1,400 × 1,050 m, cell size: 5 m) nested into a big, coarse grid (6,000 × 2,500 m, cell size: 50 m), which were linked using domain decomposition (supporting information Figure S2). Vegetation growth was only simulated on the fine grid, of which the center part (600 × 1,000 m) is used for analyses due to lateral boundary effects ( Figure S2). The morphological settings of the idealized model domain represent tidal flat systems typical for NW Europe (Roberts et al., 2000;Schwarz et al., 2014;Temmerman et al., 2007). In the cross-shore direction, the tidal flat slope decreases linearly toward the open boundary with 0.0012 mm −1 . The open boundary on the seaward side of the coarse domain is forced by a harmonic M2-tide (period: 744 min, amplitude: 1.75 m), while the other three boundaries are closed. We define hydraulic boundary conditions of 372 days, corresponding to 744 M2-tides, as 1 year. This is multiplied at every computational timestep (0.2 min) with a morphological acceleration factor (MorFac) of 30, resulting in morphodynamic timesteps of 6 min and an overall simulation period of 30 years. Each of these morphological years therefore consists of 24 M2-tides. The MorFac was chosen based on the recommended range for morphodynamic upscaling given our choices for grid cell dimensions (Ranasinghe et al., 2011) and considers the criterion suggested by Roelvink (2006) that the chosen MorFac should not significantly alter the flow field and morphological change. Furthermore, compared to previous studies with similar set-ups (e.g., Best et al., 2018;Schwarz et al., 2018;Temmerman et al., 2007), our MorFac is conservatively small.
Our simulations were initiated on a 5 m thick layer of flat, evenly distributed noncohesive sand with a median diameter (D50) of 0.1 mm and a specific density of 2,650 kg m −3 . This type of fine sand is commonly found on intertidal bars in northwest Europe prior to vegetation colonization (Schwarz et al., 2016). Sediment transport is simulated using the Engelund-Hansen equation for total transport (Engelund & Hansen, 1967;Lokhorst et al., 2018): whereby the median grain size diameter (D50) and the calibration coefficient (α) are to be specified beforehand (D50 = 0.1 mm and α = 1), while the flow magnitude (q), the Chézy friction coefficient (C), and the relative density () are cell-specific and updated continuously. Sediment availability from outside the model domain is restricted by means of an equilibrium sediment concertation at the open boundary. This set-up does not induce channel formation in case of an unvegetated marsh (see Figure S3), allowing to isolate the effects of the vegetation species.

Ecological Model
The hydro-morphological modeling was coupled to an adapted version of the ecological model by Oorschot et al. (2015) that allows to study the effect of physical plant properties as well as spatial and temporal variability in plant growth on sediment dynamics and morphological development (Brückner et al., 2019). The ecological model is coupled every M2-tide (leading to 15.5 morphological days when applying the MorFac), referred to as one ecological time step (ETS). This means that constant forcing is applied in each ETS, resulting in a stable and representative factor for imposing seasonality in vegetation, as such that 24 ETS represents 1 ecological year, and two ETS represent 1 ecological month. At every coupling moment (after each ETS), vegetation is updated based on the results of the hydro-and morphodynamic computation 10.1029/2019WR025942

Water Resources Research
and hydraulic resistance caused by the updated vegetation is altered and applied to the subsequent computations in Delft3D (Figure 1a). Simulations were run for 20 morphological years.
The effect of vegetation is incorporated in Delft3D by means of adding an extra drag term representing flow resistance in the momentum equation (-λ 2 u 2 ) based on the hydraulic resistance ( Figure 1a; Vegetation growth) calculated by the equation adapted from Baptist et al. (2007): whereby C is the net roughness, C b the real bed roughness, C d the drag coefficient, n the vegetation density (no. of stems per square meter multiplied by stem diameter), h v the height of the vegetation (m), h the water depth (m), κ the Von Karman constant (0.41), and g the gravitational force (9.81 m/s 2 ) (Hydraulics, 2006). This approach is a simplification of the way plants interact with sediment in nature where, for instance, the sediment resuspension threshold is not actively, but indirectly modeled through reduced flow velocities generated by friction of vegetation. In addition, Delft3D allows multiple vegetation types to be assigned to one cell by means of composite vegetation function (Hydraulics, 2006). This function sums the contributions in hydraulic resistance of different species scaled by their relative area cover (0-1); thus, different vegetation types occupy different fractions of a cell.
Plant colonization assumes seedlings are distributed by the tide. Thus, colonization in the model only occurs in cells of the intertidal area that have been submerged during the previous tide but are dry during low water (Figure 1a; Colonization). The actual number of cells that will be colonized depends on the species-specific probability of random establishment. Based on this value, a selection of possible cells will be colonized by the initial fraction of the respective species. The initial fraction describes the initial size of a tussock (as fraction of a grid cell) and varies among species (Table 1). If a cell is already populated by an earlier generation or other species the initial fraction reduces to the maximum available space in the cell up to a fraction of 1 (Oorschot et al., 2015). Once vegetation has colonized, spatial variability in growth is modeled through linear dose-response relationships related to physical stresses. The model incorporates mortality due to velocity and inundation (Figures 1a and S4). Mortality is modeled as a reduction in the vegetated area fraction (0-1) in a cell in response to an exerted pressure. Plants that die disappear from the model. This is based on flume experiments  where plants toppled over due to high velocities and subsequently exerted no more influence on the flow field. Species-specific sensitivity to different stresses was assigned by means of a threshold value and a slope which determines the dose-response curve. Moreover, stress tolerance was varied in two life-stages, juveniles and adults. Juveniles such as seedlings are more susceptible to stress than adult plants (Table 1). Thus, we incorporate three vegetation species in our model through different physical parameters (equation (2)) and separate stress tolerance of juvenile and adult plants (further referred to as spatial variation in growth). Furthermore, the model includes seasonal variation in growth (Brückner et al., 2019). Seasonal variation is defined from literature by species-specific colonization and growth periods and reduced plant height during winter (Figure 1b)

Species Parameters
Differences between the three considered species are incorporated in the model through differences in temporal variation in growth (Figure 1b), spatial variation in growth, and physical plant properties (Table 1). To reproduce natural development, changes in plant properties throughout ontogenesis were incorporated by assigning different properties for seedlings (first year) and mature plants (years 2-20). Since Salicornia procumbens is an annual species it has only one life stage (of one year) while Puccinellia maritima and Spartina anglica were assigned a second life stage from year 2 onward. Variability in spatial distribution of the salt marsh species are related to different degrees of susceptibility to hydrodynamic stresses: high velocities and inundation.

Water Resources Research
The physical plant properties in Table 1 are based on available literature (indicated in table) or if absent, expert judgement. Values for the chance of occurrence are partly based on observed occurrences of 5% of Puccinellia patches (Langlois et al., 2003) and the chance of plant establishment of 0.01 (Temmerman et al., 2007) and approximately 0.04 (Schwarz et al., 2014) used in modeling of Spartina. The chosen values are slightly higher because our study only considers establishment and potential connection of small patches (by adding area fractions), whereas the mentioned studies also incorporate lateral diffusion. Lateral diffusion  (Alkemade et al., 1994;Gray et al., 1991;Loebl et al., 2006), Salicornia procumbens (Davy et al., 2001;Jefferies et al., 1981), and Puccinellia maritima (Gray & Scott, 1977).

10.1029/2019WR025942
Water Resources Research was omitted in this study focusing on bio-morphodynamic feedbacks in a system dominated by seedling recruitment. We test the validity of our approach to represent scale-dependent feedbacks by determining spatial tussock distribution (random, clustered, or regular) over time for the single-species Spartina scenario, by calculating the Ripley's K (Ripley, 1977) (see Text S1). The difference between Salicornia and the other two species aims to reflect their contrasting methods of reproduction (high vs. low dispersal rates). The drag coefficients are calculated based on relations provided by Nepf (1999).

Scenarios
To investigate the impact of species assemblages on salt marsh channel emergence we simulated a variety of scenarios, classified as S1: Single species, S2: Multi species, and S3: Species shifts from one species to another ( Table 2). The S1-scenarios allow to compare the single-species effects on topography, channel characteristics, and tidal asymmetry. The S2-scenarios comprise of all possible 2-species combinations (Multi3-5) and two 3-species runs (Multi1-2). The first 3-species run uses species-specific traits (Multi1), whereas at the second 3-species run uses species-specific traits except all species have the same susceptibility to flooding (Inundation threshold, slope), which was set equal to the Spartina case (for details, see Table 1). This was done to investigate the importance of stress tolerance against inundation compared to other indirect effects that determine species colonization. Finally, two runs were performed whereby a

10.1029/2019WR025942
Water Resources Research species shift was simulated (S3). The first scenario represents a shift from Puccinellia to Spartina and the second scenario from Spartina to Salicornia. The shifts were modeled by using modeled topography of the first species after 30 years, as initial topography for the second species running for 20 more years.

Channel Extraction
Channels were extracted using a simple method relying on the difference between the initial topography and the final topography. This simple approach suffices since in our simulations only the presence of vegetation initiates channel development. Channels were assumed to be present where this difference is significant (standard deviation exceeds 0.02 m). For this extent of the domain (Ln) channels were extracted based on a negative bed level change (Z < −0.03 m) and a sufficient number of neighboring cells (to ignore isolated depressions). Based on the extracted channel networks, drainage density was calculated as the combined length of all channels over the contributing area (Marani et al., 2003): whereby δx and δy represent the grid cell dimensions and ΔY (domain width) and L n are given in number of cells. In addition, channel depth was calculated as the first percentile and bank height as the 99th percentile of bed level change.

Results
The results will be analyzed per scenario group (S1-3), while figures show the most relevant scenarios together for direct comparison. For a complete set of figures, the reader is referred to the supporting information.

Single Species (S1)
Our simulations show that the distribution of vegetation across the marsh is highly species-dependent. After 20 years of simulation, Spartina covers about 30% of the domain, while Puccinellia covers 12% and Salicornia only 1% (Figure 2a). Salicornia and Spartina grow close to the seaward boundary of the domain and are relatively evenly distributed over the range of elevations, while Puccinellia is more restricted to the upper part of the domain (Figures 2d and 2e). A slight seaward shift of the vegetation edge can be observed for both Spartina and Puccinellia (Figures 2d, 3a, and 3d), which is absent in the case of Salicornia.
The topography after 20 years (relative to the initial topography) for Salicornia does not display significant changes (Figure 3c), while channel development is apparent in the two other vegetation scenarios (Figures 3a and 3d). In the case of Spartina, a very dense and complex channel network has developed landward from x = 1,750, which coincides with the seaward boundary of the vegetation area (solid line Figure 3 a). Channels seaward of this point have lower density and complexity and more gradual levees than within the marsh (Figures 4a, 4c, and 4e). The same pattern is visible for Puccinellia (Figure 3d), yet the extent and location (channels starting around x = 2,100) of the channel network differs (Figures 4b, 4d, 4f, and S8). In addition, Spartina develops channels with a lower width-to-depth ratio (narrower) than Puccinellia (Figures 4 and 5a). Comparing total amounts of erosion and accretion per year, the maximum amount of erosion is higher for Spartina (850 m 3 ) than for Puccinellia and Salicornia (400 m 3 ). In the case of the latter two, maximal erosion occurred in year 1, while for Spartina the amount of erosion increases from year 2 onward ( Figure S9). Salicornia has a relatively constant annual erosion after 7 years have passed. Although in all scenarios the net erosion increases over the years, the Salicornia-scenario remains dominated by sediment accretion.
Furthermore, the temporal evolution of the channel network is quantified by the variation in channel depth, levee height, drainage density, and network extent ( Figure 6). All parameters increase over time following vegetation colonization and exhibit a significant difference between Spartina (higher values) and Puccinellia, although the two species initially (until year 2) show similar bed level change (i.e., levee height and channel depth) (Figure 6a). The onset of channel formation varies between the two species: year 4 in the case of Spartina and year 7 in the case of Puccinellia. For both species channel depths exceed levee heights. While Spartina and Puccinellia both cause a decline in width-to-depth ratio over time, this variable stabilizes earlier and at a higher value for Puccinellia (Figure 5a).
The presence of vegetation also influences the tidal asymmetry in domain. While the initial situation is flood dominated (Figure 7), vegetation colonization and channel emergence with Spartina and Puccinellia makes

Species Assemblages (S2)
The general development of plant cover over time of the multispecies runs (Multi1-5) is mainly determined by the presence of Spartina. If Spartina is present plant development reaches a maximum cover of about 30% (Multi1-4), yet if absent the maximum cover drops to 13% (Multi5) (Figures 2b and 2c). In the scenarios Figure 6. Temporal evolution of (a) channel depth (first percentile of bed level change), levee height (99th percentile of bed level change), drainage density (length of channels divided by surface area), and (b) channel extent (distance between upper and lower boundary of network) for S1-and subsequent S3-scenarios. The S1 Salicornia-scenario is not included since no channels developed.

Water Resources Research
where both Spartina and Puccinellia are present, the increase in cover during the first 8 years is slightly faster (Multi1,2,4) than in the single-species Spartina S1-scenario ( Figure 2).
The multispecies scenarios (S2) show a zonation in species and biomass in the main flow-direction (x-direction) for all scenarios. Since the observed plant distributions between the 2-and 3-species runs are similar, we explain our model findings exemplified on the 3-species run (Multi1). At this scenario we observe an increase in vegetation cover from seaward (x = 1,500) to the landward side (x = 2,500) ( Figure 8b). Focusing on presence/absence, Salicornia is distributed equally over a wide range of elevations, whereas Spartina and especially Puccinellia are restricted to higher elevations, visible by their convex hypsometric curves (solid lines in Figure 8c). However, when considering the area-fraction per cell occupied by each species, we find a concave hypsometric curve for Salicornia, suggesting the majority of its biomass is located at elevations between 0.2 and 0.5 m above mean sea level (0 m).
The evolution of channel networks for the S2-scenarios is strongly dependent on the combination of species present. Scenarios including Spartina and Puccinellia (Multi1,2,4) result in the largest extent of channelized area, whereas the scenario including Puccinellia without Spartina (Multi5) results in the lowest channelized extent among the multiple species runs (Figures S5c, S8b, and S8c). The presence of another species alongside Spartina reduces the overall number of channels compared to the single species runs (Multi3,4) ( Figures  S5 and S8c), as well as increases the final width-to-depth ratio (Figure 5c). The influence of multiple plant species (Multi1,2) on tidal asymmetry shows increased ebb-dominance comparable to the Spartina and Puccinellia single-species runs (Figure 7).
We tested the sensitivity of plant zonation, channel development, and tidal asymmetry in respect to different groups of model parameters. We tested the importance of spatial growth strategies (governed through inundation stress) by simulating a scenario where species only differ in their physical properties and temporal growth strategies (Multi2) compared to Multi1 where all properties were species dependent (Table 1). A

Water Resources Research
comparison between these two scenarios reveals that the total plant cover in the Multi2-scenario is about 5% higher than in the Multi1-scenario (Figure 2b). The plant distribution across the intertidal gradient is comparable between Spartina and Puccinellia at Multi2, whereas Multi1 shows a clear zonation between these two species, (Puccinellia/Spartina at high/low elevations, respectively) ( Figure 8). However, the dominance of Salicornia at low elevations remains in both cases (Figures 8c and 8d). Regarding channel development both scenarios show pronounced development of channel networks ( Figure S8b). In the Multi1-scenario, channels are more gradual/shallower with a higher width-to-depth ratio (Figure 5b). This is particularly visible at the mudflat-salt marsh boundary (Figures S5 and S6a). Regarding tidal asymmetry, Multi2 exhibits higher velocities closer to low water during ebb tide compared to Multi1 (Figures 7d and 7e).

Species Shifts (S3)
The pattern in salt marsh colonization after the species shifts is very similar to the S1-scenarios of the same species. The shift from Puccinellia to Spartina (Puc2Spar) resulted in a slightly lower vegetation cover than the S1 Spartina-scenario, while the shift from Spartina to Salicornia (Spar2Sal) resulted in a higher total cover than in the S1 Salicornia-scenario (Figure 2a).
Both species shifts resulted in a vast increase in erosion. In the case of Puc2Spar increased erosion manifested through a seaward extension of the channel network (Figures 3e and 4d). After 20 years of colonization by Spartina, channels created by Puccinellia still remained and got connected to the new marsh edge. This happened through three relatively wide and linear parallel channels (Figure 3e), which caused a rapid increase in average width-to-depth ratio across the domain (Figure 5a). Additionally, new channels formed at the vegetation edge of Spartina (Figures 3e and 4d).
In contrast, the Spar2Sal shift did not result in a significant change in channel shape or extent, but in further erosion of the entire domain (Figure 3b). This increase in erosion becomes especially clear when comparing cross-sections from before and after the shift (Figures 4a, 4c, and 4e). The Spar2Sal-scenario has resulted in deeper and wider channels as well as erosion of levees. This had little effect on channel geometry as the width-to-depth ratio only slightly reduced (Figure 5a).

Water Resources Research
Both species shift scenarios result in a rapid change in drainage density, channel depth, and a reduction of levee height just after the shift (Figure 6). In the case of Spar2Sal, the variables remain at this level for the rest of the 20 years. In the case of Puc2Spar, the levee height and channel depth increase from about 8 years after the shift onward. The drainage density reduces significantly in the third year after the shift, although it is still more than before the shift. Five years after the shift, the drainage density slightly increases for the following years.
In addition, the species shifts alter the tidal asymmetry. The stage-velocity curves in Figure 7g show a relative increase in flood-dominance in the Puc2Spar-scenario compared to the single-species Spartina and Puccinelia cases (Figures 7a and 7c). On the other hand, the Spar2Sal-scenario results in a transition from an ebb-dominant tide at the end of Spartina-dominance ( Figure 7a) to a symmetric situation (Figure 7f) that is more similar to the initial conditions.

Discussion
Our results show that channel formation and sediment redistribution of developing tidal salt marshes can be strongly influenced by the colonizing plant species. Moreover, a comparison between different species traits reveals that aside from physical plant properties (e.g., stem height, -density), spatial growth (related to resilience against inundation or currents) and temporal growth properties (i.e., phenology) have important implications on geomorphology. This is showcased by the different degrees of channel development observed between Spartina, Puccinellia, Salicornia, and combinations thereof.

Single Species (S1)
The single species scenarios showed that colonization by perennial grasses (Spartina and Puccinellia) at first did not affect the flood-dominant sediment importing system (Figures 7 and S9). However, in time as colonization by different vegetation species progressed and channels emerged, a shift toward a more (Spartina) or less (Puccinellia) ebb dominated system was observed, characterized by sediment export (Spartina, Puccinellia) (Figures 7 and S9). This was also documented in the morphologic development, exhibiting a faster increase in channel depth compared to levee height over time ( Figure 6).
To untangle the influence of growth behavior and physical species traits on marsh topography a set of additional simulations was done, where all parameters were set to a reference value except for the traits of interest. In agreement with Bouma et al. (2013), our tests indicate that the small aboveground and belowground appearance of Salicornia results in a negligible effect on topography. However, contrary to Bouma et al., who focuses on the flow field, our analysis showed no significant difference in morphological development between Spartina and Puccinellia based on physical plant properties alone. This could be linked to the resulting emerged-submerged ratio based on a M2-tide in our study compared to only two tested water levels in the flume experiment of Bouma et al. (2013). Additionally, by representing vegetation by the Baptist equation we neglect the reconfiguration of flexible Puccinellia which could be improved by alternative vegetation predictors (e.g., Cheng, 2011;Tinoco et al., 2015).
The different morphologic development between Spartina and Puccinellia in our simulations was mainly caused by different initial vegetation fractions, representing typical tussock sizes and their stem densities, both agreeing with previous studies (Bouma et al., 2009;Bouma et al., 2013;Temmerman et al., 2007). In addition, the stress tolerance of the species to inundation has a significant impact on the marsh extent and channel development. Where velocities are higher closer to the sea, the presence of a vegetation patch results in more erosion/sedimentation. This potentially caused the increased ebb-dominance of Spartina (wider) with maximum velocities at low water levels and the reduced ebb-dominance of Puccinellia (narrower) with maximum velocities at higher water levels (Figure 7). In addition, the shift from a floodto an ebb-dominated system co-occurred with the formation of wide mudflat channels which resulted in an increase in width-to-depth ratio in the Spartina-scenario (Figure 5a).
Tidal asymmetry plays an important role in residual sediment transport and the morphology of a marsh (Moore et al., 2009). Flood dominance is associated with increased sedimentation rates, infilling of the marsh and thus increased coastal protection and ebb dominance is linked to exceeding erosion of the marsh and loss of protective capacity (Fagherazzi et al., 2004;Friedrichs & Perry, 2001;Lokhorst et al., 2018). We show that different plant species cause different development of residual sediment transport More specifically,

10.1029/2019WR025942
Water Resources Research erosion exceeds accretion in both the Spartina and Puccinellia scenarios ( Figure S9), whereas accretion dominates in the flood dominated Salicornia-scenario. However, since our simulations for the sake of simplicity only use one noncohesive sediment class and restricted inflow of sediment into the domain (assumed equilibrium at the boundary), addition of fines could potentially lead to an importing system as reported in previous systems (e.g., D'Alpaos et al., 2007).
Although this study is based on an idealized domain, the vegetation growth behavior and simulated channel development are characteristic of natural systems, which typically show initial fast increase in channel density during colonization, followed by reduced network extension in the maturing marsh ( Figure 6; Allen, 2000). This is shown by simulated average rates of channel development (Spartina: 23 m/year, Puccinellia: 10 m/year) which are comparable to field observations in the Venice Lagoon (mean rate of 11 m/year; D' Alpaos et al., 2007). Moreover, van Wesenbeeck et al. (2008) found a close relation between tussock size and channel depth that generally agrees with our result of the large variation in topographic change related to (initial) fraction. Average drainage densities as calculated in this study are of the same order as previous observations in the field (Figure 6). For instance, salt marshes in Venice lagoon, Italy, range between 0.02 and 0.025 m −1 (Marani et al., 2003); Barnstable salt marshes, United States, have find an average drainage density of 0.01 m −1 (Kearney & Fagherazzi, 2016); and Salicornia-dominated San Francisco Bay salt marshes, United States, show a drainage density of 0.042 m −1 (Sanderson et al., 2000). While the simulated width-depth ratios are high (shallow channels) compared to observations in established marshes (<10, Venice Lagoon; Marani et al., 2002), they are more comparable to values based on tidal flats channels in the Dutch Wadden Sea (100-200; Marciano et al., 2005). Regarding species effects, comparison between a Spartina (Walsoorden, Western Scheldt, NL) and Salicornia marsh (Hooge Platen, Western Scheldt, NL) shows that the latter species has a weaker influence on channel formation with less deep and less pronounced channels (Schwarz et al., 2018). Likewise, Marker (1967) already observed a significant difference in topographic change in the Dee estuary between patches of Spartina (more), Puccinellia (less), and sparsely distributed Salicornia (neglectable).
However, we do not only observe comparable network properties, but also changes in the marsh edge. Since salt marsh establishment and channel development altered tidal asymmetry, it also led to redistribution in inundation times across the salt marsh. Since inundation time is one of the factors determining salt marsh's seaward location, we could observe a seaward shift of 200 m for the Puccinellia and 500 m for the Spartina marsh edge over time (Figure 2d). We expect the amount of change in marsh edge to be strongly linked to channelization of the system which offers important new opportunities to look at salt marsh development in the field.

Species Assemblages (S2)
The combination of different plant species leads to distinct zones in vegetation distribution. In the Multi1-scenario, Puccinellia is restricted to higher elevations, whereas Salicornia dominates the lowest elevations and Spartina the intermediate zone (Figures 8a and 8c). This is in line with common marsh zonation reported in literature (e.g., Beeftink, 1985;Hughes, 2004;Scholten & Rozema, 1990) and the differences in species-specific sensitivity to flooding and uprooting. Although, in our case, Salicornia was often suppressed by other species due to its low initial fraction (0.2), this is not unrealistic. Proffitt et al. (2005) observed similar situations whereby the species was suppressed by Spartina, due to reduced light availability.
The effect of species assemblages on channel development is linked to biomass and location across the tidal elevation gradient of the contributing species. For instance, during our simulations we observed a dominating effect of Spartina on channel development related to its high biomass and ability to grow low in the intertidal. However, it was striking to notice that while some species combinations (i.e., Puccinellia, Spartina, Multi4) increase the number of channels compared to the Spartina S1-scenario, others reduce it (i.e., Spartina, Salicornia Multi3) ( Figure S8) and all two-species combinations resulted in increased width-todepth ratios (Figure 5c). This suggests that competition for space modeled through temporal and spatial growth properties plays a significant role. Out of the three species, Spartina is the last to colonize, resulting in more restricted colonization due to space limitation. Consequently, cells are already filled with species possessing less hydraulic resistance (i.e., Salicornia) resulting in less channel erosion or more hydraulic resistance (i.e., Puccinellia, especially at high elevations due to less spatial mortality) resulting in more channel erosion (Figures S5 and S8c).
A sensitivity test on the importance of spatial growth properties compared to temporal strategies and physical plant properties (Multi1, Multi2) on channel development revealed that equal distribution of all three species over the tidal gradient results in less channels at higher elevations than species-specific zonation (Figures 8b, S5b, and S8b). This underlies the importance of both spatiotemporal and physical plant properties in modeling salt marsh ecosystems, and shows that although ecological interactions (e.g., competition) are not considered in the model, indirect hydromorphodynamic-plant interactions do result in realistic plant assemblages (Multi1). But it also opens new questions in how vegetation zonation might be able to optimize drainage efficiency of intertidal marshes, which was previously only linked to vegetation presence-absence (Kearney & Fagherazzi, 2016).

Species Shifts (S3)
We investigate the effect of species-shift scenarios probable to occur over climate change on channel development. Marshes are believed to naturally respond to an increase in sea level by means of increased sedimentation (Kirwan & Megonigal, 2013;Morris et al., 2002). However, on the long term, sea level rise is accompanied with transgression, meaning that upper marsh species are replaced by lower marsh or pioneer species that are less susceptible to inundation (Puc2Spar). Moreover, as shown in the Po-delta in Italy, rising temperatures and increased drought benefits Salicornia over Spartina pioneers (Spar2Sal) (Strain et al., 2017). In this study we aimed to isolate effects of such climate-induced species shifts rather than look at more direct climate effects (e.g., sea level rise).
Shifts between colonizing plant species had an important effect on the resulting channel networks and sediment distribution. The shift from Puccinellia to Spartina (Puc2Spar) resulted in rapid extension of the channel network, which developed into wider and deeper parallel channels, because Spartina is less sensitive for inundation and grows much more seaward (Figures 3 and 4). This type of channels allows increased ebbdominance and connectivity between the marsh and the adjacent open water body reducing the protection of the hinterland (Figure 7) (Leonardi et al., 2018). The shift from Spartina to Salicornia (Spar2Sal) resulted in wider channels and lowered levees (Figure 4) in support of the theory that vegetation has a stabilizing effect on channel banks. The Puc2Spar shift results in a reduced width-depth ratio of the channels (D'Alpaos et al., 2006;Schwarz et al., 2014), while the transition from a densely populated marsh to an almost bare mudflat causes channel banks to erode. Consequently, a species shift from a perennial grass to an annual halophyte (or almost bare mudflat) would be disadvantageous for the protective capacity of the marsh. This has also been pointed out by Strain et al. (2017) who state that an increasing dominance of Salicornia (veneta) would eventually result in a marsh with reduced resilience, and a reduced capacity to respond to sea level rise.
Both species-shift scenarios show that initial channel configuration has an important consequence for the final network. This is showcased by the fact the neither the Puc2Spar nor the Spar2Sal are comparable to any of the other single species scenarios (Figures 3 and S8). We implemented the species shifts as an abrupt shift which might be more gradual in nature. The increase in drainage density that companied with the species shift is an artefact of the way drainage density is calculated and should be omitted (equation (3)). However, we could show the importance of climate-induced species shifts on marsh morphology, which need further field validation in order to predict all their implications on coastal management and protection.

Conclusions
The aim of this study was to investigate the impact of different plant species on salt marsh and channel development by means of a bio-geomorphic model study. We specifically investigated the effect of single species (S1), multiple species (S2), and species-shifts (S3).
The single species scenarios (S1) indicate the influence of vegetation on the marsh topography to be highly species dependent, which is related to spatiotemporal variations in vegetation fraction and physical plant properties. Spartina and Puccinellia induce significant channel formation, while Salicornia did not induce 10.1029/2019WR025942 Water Resources Research channel formation. The multispecies scenarios (S2) show that the species in the assemblage influence channel formation depending on their biomass, and that different combinations facilitate or hamper channel development compared to single species scenarios. The species shift-runs (S3) show that channel network characteristics and sediment transport are very different from other runs underlining the importance of climate-induced species shift on morphology. The results reveal that the medium-scale morphological development of tidal flats strongly depends on salt marsh species and composition. Climate-induced changes in salt marsh composition will lead to strong effects on the channel network and drainage capacity on the marsh.