Distinct patterns of bank erosion in a navigable regulated river

Distinct bankline patterns appeared after the removal of protection works along a navigable reach of the Meuse River. A series of oblique embayments now dominate the riverine landscape after ten years of bank erosion, but their location and asymmetry cannot be explained yet. This work analyses and integrates field measurements of flow, ship waves, bank composition, bed topography and historical maps to explain the observed patterns along two reaches of the river. An extraordinary low‐water‐level event generated by a ship accident provided the unique opportunity to also analyse the subaqueous bank topography.


Introduction
Since the EU Water Framework Directive (WFD, 2000) several measures have been taken to improve the water quality and habitat diversity of European rivers (Kallis and Butler, 2001;Pearce, 2013;RESTORE, 2018). Particularly for riverbanks, bioengineering protections and complete protection removals have already been implemented in several countries, such as Austria (Liedermann et al., 2014), Germany (Schmitt et al., 2018) and France (ONEMA, 2018). In the Netherlands, the Meuse River constitutes a unique case where 40 km of banks had their protection works removed between 2008 and 2016. The river functions as an important waterway with water levels regulated by weirs. The implementation of this large-scale bank re-naturalisation project comprises different phases. The first one, finalized in 2018, serves as a reference for the next phases until 2027, when the restoration works will conclude. In this context, the monitoring and understanding of the morphological evolution and ecological development of the restored banks is crucial to gain knowledge to improve the impacts of further operations.
Ten years after the first reaches were exposed to bank erosion, the new channel margin appears non-uniform at several locations, where it is often characterized by the presence of oblique embayments. In particular, two reaches near the city of Gennep present either deep embayments (Figure 1) or uniform bank alignment. By analysing sequential airborne laser scanning, Duró et al. (2018a) showed that the sediment yield varied greatly in these reaches due to wide-ranging erosion rates, and that the erosion patterns formed above a terrace or bench, consisting of a flat surface with elevation close to the minimum water level. Yet, the factors determining the terrace and observed patterns of erosion remain unknown.
Heterogeneous bankline patterns occur at different spatial scales. At the scale of the channel depth, irregularities appear due to mass failures and slumping, leaving irregular sediment deposits at the bank toe (Darby et al., 2010), or due to roots that create and support protrusions (Rutherfurd and Grove, 2004). These features are typically associated with variations in bank roughness (Kean andSmith, 2006, Darby et al., 2010). At the bend or bar scale, banklines are curved, sometimes sharply, associated to bar and bend dynamics (Thorne et al., 1993;Klösch et al., 2015) and secondary flow (Kleinhans et al., 2009;Ottevanger et al., 2012). Bank irregularities at an intermediate scale are not common, and none of the mentioned physical processes can explain their formation. Finally, embayments of similar size to those in the Meuse River appeared in the Mekong River (Hackney et al., 2015), but the size of the latter is substantially larger than the former.
The causes of the erosion patterns observed along the Meuse River are not easily discernible. Tree roots typically delay erosion by increasing the resistance against bank instability (Pollen and Simon, 2005;Wiel and Darby, 2007) and by reducing effective shear stresses over the soil when protruding in the flow (Khanal and Fox, 2017). However, not all embayments are delimited by trees. Ship-induced waves exert oblique loads onto banks (CIRIA, 2007), but the directions of embayment evolution are not consistent with ship movements. Composition and stratification of riverbanks do affect erosion rates and bank failure (Thorne and Tovey, 1981;, whose variabilities are also found along single meander bends (Konsoer et al., 2016b). Yet, no evidence of longitudinal changes of soil compositions have been reported to create oblique embayments, to the best of our knowledge. The formation of a terrace abutting banks in rivers used as waterways has been reported (e.g., Liedermann et al., 2014;Murray-Darling Basin Authority, 2017), but the contribution of ship waves to erosion is usually unclear due to other factors that may play a role, such as floods , sapping (Hagerty et al., 1995;, or wind waves . This is partially due to the difficulties to quantify small amounts of erosion produced by single ship passings (Bauer et al., 2002), for which prolonged monitoring of erosion is necessary for accurate yet cumulative measurements.
The objectives of this study are to unravel the drivers of the bank erosion patterns that have appeared along the Meuse River and to provide insights on their contribution in regulated navigable rivers. We hypothesise that ship waves, floods, trees and bank composition have different roles on the formation of the bankline patterns, and thus analyse each factor to disentangle the dominant ones. The method adopted here is based on the analysis of measured data, field observations and historical maps. Notably, an extraordinary low-water-level event generated by a ship accident allowed us to analyse also the subaqueous bank topography in detail.

Study sites
The Meuse River is characterized by a pluvial regime with peak flows in winter/spring reaching 3100 m 3 /s, and discharges during summer/autumn as low as 40 m 3 /s (Descy et al., 2009). In the past, the Meuse used to meander across its floodplains . The river was trained during 1940s to 1960s to facilitate navigation and increase water conveyance. Poplar trees were planted every one hundred metres along the banks to guide ships during overbank flow conditions. The main channel was then canalized with a regular cross section (120 m wide, 1:2.5 bank slope), bends were cut off, and several weirs with ship locks were constructed (e.g., Sambeek is visible in Figure 2a). The relatively recent restoration removed several of those rip-rap revetments and transverse groynes, to promote natural processes and riverine habitat diversity.
The study reaches are located near the city of Gennep (Figure 2), the Netherlands, where the current riverbed is composed by sand and gravel and the floodplains are covered by grassland and by agricultural fields further away from the main channel ( Figure 1). On both sites, the protection works were removed only along the left bank ( Figure 2). Reach A was restored in two stages, first in 2008 along 750 meters and then in 2010 over the upstream 250 m. Reach B extends over 1.2 km downstream of reach A and was restored in 2010. The bankfull depth at both sites ranges between 10 and 11 metres, with minimum regulated depths of c. 7 m.
Since the restoration works, banks retreated in both reaches at wide-ranging rates, resulting in distinct erosion patterns . Notably, the largest embayments along each reach had either consistent upstream or downstream asymmetrical orientations. Some bank stretches presented a rather uniform bankline and others show symmetrical embayments with somewhat smaller dimensions than the asymmetric ones. All stretches evolved forming either parallel banks or embayments over the years. In addition, a bench at the bank toe was present all along these reaches, showing that the bank retreat mainly occurred above the terrace level (Duró et al., 2018a(Duró et al., , 2018b.

Methods
We examined the role of factors which could produce the observed bankline retreat patterns through analyses of available data. We first studied the capability of flow currents and ship waves in producing a load distribution that could match the observed bank retreat patterns over the years. Then, we analysed the bank strength through its composition, at a local scale with deep cores and a larger scale with surface samples, and we considered the presence of trees. We examined in detail the mechanisms of upper-bank erosion over a year, and investigated past river dynamics to infer processes related to lithological characteristics. Finally, we considered the contribution of other factors such as groundwater sapping erosion, cattle and rainfall-induced erosion.
We inferred  what can have caused the bank patterns, distinguishing between initial conditions, such as inherited geology or planted trees, and mechanisms, such as rain-or ship-driven, spatially-varying erodibility. The role of currents, ship waves and bank composition were analysed through field data of reach B, which are described in detail in the following sub-sections. The channel migration in the last period preceding canalization was studied through historical maps and recent aerial photos. The influence of trees, sapping, cattle and rainfall was evaluated based on field observations along both study reaches. The data used for the analysis came from three different sources. First, the Ministry of Infrastructure (Rijkswaterstaat) provided measured and validated time series of discharge at Venlo and of water levels at Gennep and Sambeek, the river bathymetry and yearly aerial photographs. Second, we surveyed the terrace and upper-bank topographies, measured ship-induced waves, collected and analysed soil samples, performed deep cores in the floodplain, and took photographs and videos from the field sites. Third, the National Archive of the Netherlands made available historical maps where the case studies are located, which was facilitated by courtesy of Rijkswaterstaat. All elevations in this work refer to NAP (Dutch reference sea level).

Mechanisms for erosion
We cross-analysed the discharge and water level series and the timing and magnitude of average current-induced bed shear stresses at a location with active erosion and relatively large bank retreat. For that, we took the channel cross section at km 153.940 of reach B (Figure 2a) as reference, due to the relatively high erosion rates during 2017. The average bed shear stress was calculated as assuming hydraulic radius equal to average water depth (h) since width/depth > 10, with ρ=water density (kg/m 3 ), g=gravity acceleration (m/s 2 ), S 0 =average energy slope. Energy slopes and water depths were computed based on a linear interpolation between the two nearest known water stages. An acoustic Doppler velocimeter (ADV) placed at km 154.0, before the beginning of the terrace measured water levels for a period of three weeks, with a frequency of 8 Hz to capture short waves. We identified typical wave characteristics, such as period and height of both primary and secondary waves. We also observed the generation, propagation and breaking of ship waves across the channel and over the terrace in reach B during several field visits, and analysed their relation with the morphological features of the bank area. In order to compare the order of magnitude of current-and ship-induced shear stresses, we estimated wave-induced shear stresses on the bed at the start of the terrace with the linearwave theory and Jonsson's (1966) concept of wave friction factor f w , as explained next.
The maximum bed shear stress induced by secondary waves using the friction factor results in with u=velocity amplitude near the bed (m/s). The amplitude of the near-bed velocity is with ω=angular frequency=2π/T (1/s), T=wave period (s), H=wave height (m), k=wavenumber=2π/L (1/m). The wave length L (m) is obtained through the dispersion relationship The friction factor as expressed by Swart (1974) with a threshold of 0.3 when (ξ/r) < 1.59; ξ is the particle displacement equal to ξ =U 0 /(2π/T 0 ); r is the bed roughness taken as 2 cm for the Meuse to account for bottom irregularities. The upper bank topography was measured using an unmanned aerial vehicle in combination with structure from motion photogrammetry (Westoby et al., 2012;Clapuyt et al., 2016), following the methodology of Duró et al. (2018b) to measure riverbanks. We applied this methodology in nine surveys along 2017. We used a DJI Phantom 4 with 18 ground control points to georeference the model and Agisoft PhotoScan to process the imagery. The resulting digital surface model achieved a resolution of 2 cm and a root-mean square error of 3 cm by comparison with RTK GPS points.

Factors affecting bank resistance
We measured the subaqueous topography down to 2.5 m below the minimum stage and took seven surface samples along reach B ( Figure 2) on January 2017. Following the studies of Kimiaghalam et al. (2016), we used cohesion (C) as an indication of the material erodibility to make a hierarchy of erosion resistance between different areas. Then, we estimated critical shear stresses for entrainment (τ c ) through the linear relation with units of Pa for τ c and kPa for C. The digital topography of the terrace was used to identify lithological layers along the reach and measure their ridges, which were validated through observations in the field. The stratum ridges that presented varying elevations were not included in the analysis because their superficial appearance might be shaped by erosion, not corresponding to actual stratum strikes and dips. Those banklines at large embayments that have asymmetrical orientations with respect to the channel axis were also marked.
Furthermore, stratum dips were also computed by the theory of Struiksma et al. (1985) adapted by Talmon et al. (1995), similarly to van de Lageweg et al. (2014). The prediction of the dip for lateral accreting deposits in an infinitively long bend reads: where z=bed elevation, n =transverse direction in a curvilinear channel for damped conditions, θ=τ/((ρ s -ρ)gD) Shields mobility parameter, κ=von Karman's constant, C=Chezy coefficient (m 0.5 /s) calculated as 18*log(12 h/D 90 ), and R=bend radius of curvature (m). The material at the terrace surface was analysed at 7 locations, chosen to cover areas of different erosion magnitudes, from the least to the most retreated banks (see later Figure 4). Coring was done with a 15 cm-long cylinder and after removing the top 10 cm. All surface samples were subjected to direct shear tests to derive the internal friction angle and cohesion of the soil by least-squares linear regression over a range of normal loads of 17, 36 and 73 Pa and shear rates of 0.01 mm/ sec. Later, all samples were dried in the oven at 105°C and mechanically sieved and weighed to obtain the particle-size distribution. For all fractions smaller than 63 μm the granulometric curves were extended to 2 μm by hydrometer analyses to distinguish silt from clay fractions, measuring relative density changes of water as the mud settles. All test were done according to the standard: BS 1377-2:1990.
In total three deep corings were performed near locations were both minimal and large bank retreat has occurred (Figure 4). Borehole levels were measured using a RTK GPS. In total two different coring techniques were used to retrieve sediments namely i) Edelman corer (above groundwater levels) and ii) Van der Staay suction corer (below groundwater levels). Sediments were logged in the field at a 10 cm interval in terms of lithology (USDA classification) and other sedimentary characteristics (Berendsen, 1982).
The river dynamics during the last years of free migration were analysed through two historical maps that indicate positions of the main channel prior to the canalization. We used a map from the preliminary project for the Meuse River canalization (Nederlandsch-Belgische Commissie, 1912), which also indicates the extent of a large flood that occurred during 1880. We then used a second map dating from 1950 and belonging to a series of national river surveys (Netherlands Nationaal Archief, 1952). Both maps were georeferenced in ArcGIS with several landmarks that served as fixed points in time, such as churches, roads, and bridge piles. We then digitized the channel banklines at all distinguishable angle changes.

Currents
Both study sites are rather straight and submerged bars are absent in the channel, meaning that bend flows are negligible (Papanicolaou et al., 2007). The pluvial regime of the river is evident in Figure 3a through the winter/spring peak flows that stand out over the summer/autumn low discharges. Since 1940s, water levels have been strongly regulated to ensure navigability, which results in regular and extended periods of rather constant stages during low discharges. As a consequence, bed shear stresses become significant only during floods, i.e. for a relative short time during the year. Therefore, current-induced shear stresses at the bank are expected to have similar durations but lower magnitudes, as for instance, 60-80% of the bed shear stress (ASCE, 1998). Figure 3b presents the river a cross section located in the longest embayment of reach B. The left bank had the rip-rap protections removed in 2010. At that time, the bank profile had a 1:2.5 slope from the toe up to the floodplain level. Currently, the bank presents a terrace at an elevation of 7.0-8.0 m that extends over 20 m inland, before encountering a 3.5 m bank scarp. The submerged part of the bank presents a c. 3 m high notch that was dredged in 2012 for ship manoeuvring, whose remains partially rest at the toe. The right bank is still protected by rip-rap and belongs to a breakwater that divides the main channel from the harbour (Figure 2a).
The water level frequencies in the period 2009-2017 are plotted on the right axis of Figure 3b, overlaying the river cross section. For more than 70% of the time, the water level fluctuated within a range of 0.50 m, which coincides with the elevation range of the terrace (7.5-8.0 m). The terrace extends over the entire reach and defines the toe of the bank scarps. The bankline patterns visible in Figures 1 and 2 formed above this bench. Water level frequencies also show that overbank flows rarely occurred in 7 years. Since 2009, the floodplains were inundated only at the peak of the largest flood event in January 2011 (Figure 3a).

Ship waves
Navigation in the study reaches, located between Sambeek and Grave ship locks, accounts for about 35,000 passings of commercial vessels and recreational boats every year. The water level measurements indicate that primary waves have periods ranging from 25 to 65 seconds and amplitudes up to 0.45 m. Secondary waves have shorter periods, namely 1.25 to 3 seconds, with typical values of 2.25 seconds and amplitudes of 0.10 m. The maximum recorded amplitude of secondary waves was 0.45 m. According to linear wave theory, the low but frequent waves with height of 10 cm and period of 2 seconds induce a maximum bed shear stress of approximately 0.6 Pa through orbital velocities. As water stage is highly controlled (Figure 3b), the generated waves mostly impact banks at a narrow range of elevations.

Lithological succession and stratification
The topography of the terrace has a mild slope towards the lower bank and presents stratification at certain locations, exposed after the erosion of the bank material above (Figure 4). A downstream view from km 153.950 clearly shows the  (Figure 5a). Here, the sedimentary layers crop out as minor scarps and form an angle of 12 degrees with the main channel axis, similar to the orientation of the bankline next to it. At other embayments, the layer orientation is more similar to those of the banklines (see dashed lines in Figure 4c), laying approximately parallel to each other. Ridges in Figure 5a are more irregular in the background, probably due to the longer time this area was exposed to erosion. Figure 5b presents an elevation profile of the terrace across layer strikes where some of the strata in Figure 5a are recognizable. In this profile, the strata are inclined by 5 degrees with respect to the horizontal plane. The predicted dip based on Struiksma et al. (1985) ranges between 4 and 6 degrees, in agreement with the measured one. This range considers a constant water depth of 8 m, a Chezy coefficient of 59 m 0.5 /s with D 90 =0.05 m, and a variation of the Shields parameter between 0.7 and 1.8. The θ range results from either considering τ=1 Pa and D 50 =0.09 mm or τ=2 Pa and D 50 =0.07 mm, which account for possible values of shear stresses during high flows and mean particle size of soil samples (see next).
The samples along reach B present wide-ranging compositions at the terrace level. Clay contents range from 8% to 25% and silt from 23% to 75%, thus sand content is also diverse, ranging from 9% to 68% (Table I). Such heterogeneity in grainsize results in different textures according to the USDA soil classification. Samples 1, 2 and 4 are sandy loams, 3 is a loam, and 5-7 are silty loams. In addition, critical shear stresses for entrainment, which are linearly related to cohesion, follow a trend with the texture classes. Sandy loams present lower critical shear stresses than silty loams and loams.
The location of samples with the highest cohesion, and thus highest critical shear stresses, correspond to the least retreated banks. Sample 3 is located between deep embayments ( Figure 4b) and belongs to a protruding but submerged soil layer, higher than the surrounding terrace elevation. Samples 5, 6 and 7 belong to the uniform stretch of reach B (Figure 4a).
These four samples have C >14 kPa and τ c > 12 Pa. Samples 2 and 4 belong to deep embayments, and sample 1 to an area with intermediate retreat. These samples have C ≤ 13 kPa and τ c ≤ 11 Pa, corresponding to the lowest range of the sampled materials. Consequently, there is a reasonable correlation between relative erosion rates at key areas of reach B and the compositions of the different layers at the terrace level, especially regarding their cohesion.
The lithology around the embayment of Figure 4b also presents varying classes throughout the floodplain depth ( Figure 6). The upstream core located next to the oblique bankline, number 3, mostly displays silty-clay loam from the bank top down to 9.0 m. The next 2 metres are mainly composed by loam. From 7.0 m down to 5.3 m, there is a mixture of loam and sand layers. Below 5.3 m, the core shows sand with some gravel contents. Core 2, located in the direction of the embayment evolution (Figures 2d,5b), presents variations between silty loam and silty-clay loam down to 9.0 m, followed by a loam layer of 0.5 m. Then, from 8.5 m down to 6.2 m, there are several sediment layers ranging from sand to loam. Deeper than 6.2 m, the lithology is mainly sand with some traces of gravel. Core 1, located near tree 4 and the least retreated bank area, has a mixture of silty loam and silty-clay loam at the top bank down to 8.5 m. The following 2 metres display loam and sandy loam layers, the former having traces of organic matter. From 6.5 m down to 5.4 m, the main composition is loam, presenting some variations within thinner layers and traces of organic matter. Below 5.4 m, sand is ubiquitous with the exception of thin loam layers and traces of gravel.

Upper-bank erosion
Between February and November of 2017, bank erosion at km 153.940 progressed uninterruptedly (Figure 7). Between the first two surveys there was a flood event having a duration of  two weeks (Figure 3a) that raised the water level up to 9 m. During the rest of the time, the water level remained at 8 m.
The difference between the first two topographic surveys show erosion throughout almost the whole bank height. From March 15 to April 26, only toe erosion occurred. Afterwards, a mass failure happened. Subsequently, from June 8 until August 23, slump blocks and the toe progressively eroded, reaching an incipient undermining of the bank. Between August 23 and October 11, the bankline retreated and the wasted material was removed. The last survey on November 23 shows minor toe erosion.

Currents
The flow structure in the near-bank area can be complex and three-dimensional (Rhodes and Knight, 1994;Blanckaert et al., 2010), particularly at the sharp expansions of embayments. Here, the boundary layer could detach (Simpson, 1989), triggering a different type of erosion mechanism, as scour holes in river beds develop due to turbulent mixing layers (Hoffmans and Booij, 1993). Similarly, embayments could grow in the horizontal plane from initial bank irregularities due to recirculating flow (Hackney et al., 2015).
Since 2009 floods occurred over relatively short periods, with water levels below bankfull. The mechanism of vortex shedding within an embayment needs to act at different spatial scales to produce the erosion patterns, from initial stages of bankline retreat until full development of the embayments (Figure 2). Yet, this hypothesis neither explains different directions in embayment evolution nor the formation of the patterns only above the terrace level. Furthermore, flow-induced bed shear stresses at peak discharges are lower than the critical  shear stresses of all soil samples (6.4-18.5 Pa). Hence, the contribution of flow currents to embayment growth and terrace formation was likely minor.

Ship waves
The estimated bed shear stress induced by propagating secondary waves of average height is smaller than the flow-induced bed shear stress at peaks (0.6 versus 2.0 Pa), but occurs much more frequently. However, the assumption of linear wave theory leads to underestimated shear stress, because ship waves are steep and shear stresses exceed soil entrainment thresholds. The latter is indicated by the plumes of suspended sediments that originate from the banks during low flows (Figure 1). Therefore, waves explain the gradual terrace advance during low flows (Figure 7).
Moreover, the planform evolution of the restored banks presents four characteristic types of bank retreat that could be related to the load distribution exerted by ship waves. Banks retreat a) parallel to the channel centreline (e.g., between trees 4-5, Figure 2b, or trees 10-11, Figure 2e); b) creating embayments that grow without a clear trend towards neither upstream nor downstream (between trees 3-4, Figure 2b, or trees 9-10, Figure 2e); c) creating embayments that evolve asymmetrically towards either upstream or downstream (those with arrows in Figures 2b-e); d) creating embayments that evolve asymmetrically and with a certain degree of irregularity (between trees 4-5, Figure 11d). All these types are schematized in respective panels of Figure 8.
Considering that ship waves dissipate over the distance as they propagate, the longer they travel the lower energy they carry. In addition, a terrace with shallow water creates higher resistance for wave advance than deeper areas, especially if waves break. These considerations promote faster bank erosion at least retreated areas (Figure 8a) than inside the embayments (Figure 8b), where longer and shallower areas are present and waves refract and diffract lowering the specific energy. Yet, this contradicts the greater bank erosion rates observed at embayments compared to uniformly retreated stretches (see bankline evolution after two years of restoration, Figure 2b-d), which is coherent with the evidence of active erosion at the embayment presented in Figure 7 after six years of restoration.
At the terrace and in the presence bay-shaped banklines, waves refract or diffract before they break. Sharp and smooth bankline changes influence wave propagation and in turn bank retreat, which may result in non-linear interactions affecting bay formation. Once an initial perturbation grows into an embayment, asymmetric evolution (Figure 8c) can be related to primary and secondary waves. Deep primary waves create strong localized currents along sharp bankline changes at embayment extremes, acting like bores and locally increasing the load due to momentum change along their path.
In addition, secondary waves may enhance bays at extremes since here they approach banks at approximately right angles (see Figure 8c). The asymmetric evolution could thus be caused by a higher number of loaded ships in one direction (e.g., downstream) due to higher waves impacting a given bay extreme. However, this cannot explain different evolution directions of the embayments in reaches A and B, both located on the left riverbank.
The last type of bank evolution ( Figure 8d) presents a change in the rounded shape of embayments that results in sharp bankline change. This particular morphology does not follow the typical evolution of other embayments, for which it can be related to local conditions (dashed box in Figure 8d). These are elaborated in the next sections, including strong root systems (Figure 2d, tree 5) or low erodible sediment layers. To conclude, ship-induced waves are incapable of a explaining the different bankline retreat modes in a comprehensive way, despite concentrated loads at embayment extremes may promote their growth.

Trees
Root growth in fluvial environments is mostly dependent on water and oxygen availability driven by water table fluctuations, which are conditioned by the river flow regime (Rood  , 2003;Tron et al., 2015). At the case studies, all trees grew approximately 30 m above ground level and developed their roots under the same weather and water table fluctuations, which are particularly limited by the minimum regulated level in the river. Since these controlling factors for root growth were similar for all trees, then similar root size and structures are likely among them. Figure 9 shows trees 10 and 11 of reach B at two different stages of erosion.Despite the presence of extended roots in the upper bank, the response to erosion of trees was disparate. On the one hand, the surroundings of large embayments ( Figure 1) present areas with relatively low retreat, which coincide with the location of the 30 m high trees, such as trees 4 and 5 in Figure 2d. This could result from rootreinforced soils and associated chemical strengthening by extra-cellular polymeric substances. On the other hand, in other areas erosion surpassed the location of trees (e.g., trees 7 and 10 in reach A, and trees 3,6 and 7 in reach B). Also, some embayments are not confined by trees, as tree 8 is 40 m from the end of the largest embayment of reach B (Figure 2e), whose bank scarp does not present any roots. As a consequence, the presence of trees in the upper bank does not offer an exhaustive explanation for the planform patterns of erosion.

Lithological succession and stratification
The floodplain indeed presents highly heterogeneous compositions arranged into tilted strata along the reach, whose strikes are oblique to the current channel position. Figure 10 shows two channel positions before the canalization works of the 1960s to infer possible depositional processes during floodplain formation (e.g., Lewin and Ashworth, 2014). In 1912, the river main channel presented two clear meanders connected to a rather straight reach (Figure 10a). The latter bend migrated from 1912 until approximately 1950 smoothing the curvature while shifting downstream and to the right (Figures 10c-d), involving processes of erosion along the outer bank and accretion at the inner bank. The channel was then fixed with revetments and groynes, and later canalized cutting the bends off, resulting in the present channel alignment (dashed lines in Figure 10). Figure 11 compares the orientations of the inner banklines in 1912 with the layer strikes identified on the terrace (Figure 4), including also the current bankline orientation of asymmetric embayments.
The orientations of strata and the 2017 banklines at embayments match the orientations of the 1912 inner bank along the reach, except for an outlier at km 153.6, which could be the consequence of a local discontinuity in the bank resistance, forming a mild abutment. The series have an average orientation difference of 11 degrees and a 45 metre downstream shift. This correspondence suggests that the meander bend migrated to the 1912 position passing through the current channel location, depositing sediments in sequences that determined the observed strata. What is more, the channel migrated from 1912 to 1950 with a similar drift towards the downstream-right direction, also with rotation and translation of banklines from an absolute reference system. This downstream shift of the inner bankline can be interpreted as the result of scroll bar deposits Candel et al., 2018). The measured dip of the strata is  Struiksma et al. (1985), for lateral accreting deposits of meandering rivers (Gibling and Rust, 1993). Moreover, the deep cores around the embayment present lithological successions that vary in elevation and thickness, which is explained by successive scroll bar formation. Figure 12 interprets the ridge-swale morphology from the deep cores and illustrates the planform layout of former and current channels.
The depositional sequence follows the direction of the past river migration, which is oblique to the current channel position, so embayments appear as a consequence of these varying compositions along the present channel (Figure 12b). Importantly, the elevation of the point-bar deposits especially varies at the water level range (7-8 m) where primary and secondary waves normally attack the banks (Figure 12c). Above these variations, overbank deposits fill the top bank with fine sediments, dominated by the presence of loam up to 9 m and above that by silt, levelling the floodplain morphology (Nanson and Croke, 1992). These spatial variations agree with the different lithologies observed at the terrace along reach B (Table I). Therefore, bank retreat rates and patterns are controlled by the nature of the deposits, producing faster erosion at sandy deposits and lower at loams, whose spatial disposition is defined by the structure and orientation of strata.
The asymmetry of embayments is caused by the obliqueness of sedimentary layers with respect to the canalized channel. The consistent but opposite asymmetric orientations of embayments in reaches A and B correspond to scroll bar formation in different directions during meander migration (Figures 2a,10a). On the other hand, parallel bank retreat, such as the straight bankline in the downstream end of reach B, responds to a relative strong layer at the minimum regulated level laying almost parallel to the current channel position (Figure 10c). Finally, the different embayment lengths are most likely the result of varying layer thicknesses and strike angles with the main channel, which create a variation of the projected length along the channel (see Figures 12b,c).
The layer that protrudes at the upper bank near tree 4 exemplifies the control exerted by lithology and stratification on bank retreat (Figure 13a). This layer with relatively high clay content (sample 3, Table I) is aligned with the subsequent downstream bankline on the right, corresponding to the beginning of the next embayment ( Figure 4). Its oblique orientation with respect to the canalized channel defines the asymmetry of the downstream bay. Even though tree 4 locally delays erosion rates by root reinforcement, the response of bank retreat follows the location of this cohesive layer (Figure 13b), except where shear stresses are highest due to less wave dissipation and high flow currents. The upstream bay results from the erosion of the scroll-bar ridge composed by sand at the regulated water level (sample 2 and Figure 12), whose development rate reduces at the encounter with the downstream scroll-bar swale, filled with loam at the controlled water stage.
The disparate erosion at the location of trees is then explained by the erodibility of the layers at the upper-bank toe, and the relative position of the trees across them. Figure 13b schematizes the horizontal extent along which a strong layer delays erosion rates in front of a tree, called 'survival range' to refer to its short-term fate. This range depends on the layer thickness, tilt and position with respect to the tree, but also on the elevation range of wave action (between dotted lines) propagating at controlled stages (upper limit of green area). Those trees over a weak stratum are dislodged after few years, while those partially covered by a strong one (e.g., tree 3 of reach B) are later outflanked.

Other factors
No evidence of sapping was observed along the hydrological year in the study reaches, neither in the form of regular cavities at bank scarps nor through deposits of dislocated particles (Hagerty, 1991). Cattle grazing on floodplains may widen small streams by progressively breaking banks down (Trimble, 1994), but in the study sites banks heights reach 3.5 m and cow pathways were always observed far enough from the edge to avoid mass failures. Rainfall events can increase the soil specific weight and induce bank mass instability (Simon et al., 2000), but only a few local and isolated failures were observed at the beginning of the rainy season.

Ship-induced waves and terrace formation
The relative contribution of waves on eroding banks depends on vessel frequencies and characteristics and on the natural forces also acting on the banks. The case studies present highly regulated water levels that allow ship waves to attack banks at a narrow range of elevations. These conditions are sufficient to Figure 12. Scroll-bar stratigraphy interpreted from deep cores (schemes adapted from Gibling and Rust (1993) and  develop a terrace across riverbanks during low flows, despite the erosion produced during floods. Furthermore, the fact that the water does not surpass a minimum stage impedes some mechanisms of erosion to act on the lower bank, such as drawdown (Simon et al., 2000) and subaerial processes (Wynn et al., 2008), which further increases the differential retreat between the lower and upper bank. Detailed monitoring of the upper-bank erosion showed that floods were not necessary for the basal clean-out of failed material , and ship waves acted disaggregating and removing slump blocks, as observed in the Ohio River (Hagerty et al., 1995). The role of ship waves extend over the whole erosion cycle, including undermining and destabilization of the bank top. In the most erodible stretches, the terrace presents mild slopes of c. 1:25, unlike situations with more frequent water level changes that display steeper bank profiles . Therefore, highly regulated rivers with high ship traffic set the conditions to form a mild extended terrace across banks, which reduces erosion rates over time due to wave energy dissipation in shallow waters, but continues developing after 8 years of protection removal.

Effects of tree roots on bank erosion
Uniformly vegetated banks can significantly reduce near-bank velocities along entire river bends (Konsoer et al., 2016a). At a smaller scale, Pizzuto et al. (2010) suggested that the joint effect of nearby trees on flow detachment could reduce erosion rates similarly to small-scale roughness (Kean and Smith, 2006), based on short-term observations of sequential abutments on sandy-loams. Rutherfurd and Grove (2004) demonstrated that isolated trees locally delay erosion rates in sandyloams, which is shown by root-plate abutments in the bankline, but their effect is negligible for the migration rate of the meander bend, in agreement with our results.
The combination of large woody roots and different substrates led to disparate erosion resistances. The effect of trees located on highly erodible layers in the study reaches (e.g. tree 7 over sandy loam, Figure 4) appeared negligible at yearly intervals of bank retreat, even at a local scale. On the other hand, trees with cohesive substrates presented high erosion resistance. Vannoppen et al. (2017) showed that the additional resistance against concentrated flow erosion of soils with thick taproots increases as the cohesion increases, but decreases with increasing sand content. Hence, the mechanism through which woody roots reduce bank erosion rates is more effective with less erodible substrates.
This positive feedback between the presence of vegetation and cohesive soils during the entrainment phase of the erosion cycle (Thorne and Tovey, 1981) can explain the divergent fates of trees in the study cases. Cohesive soils permeated by deep roots hold steeper banks (Thorne, 1990) by increasing the resistance against instability (Pollen-Bankhead and Simon, 2009). This delays the failure mechanism compared to the case without root reinforcement, allowing for a longer phase of entrainment within the erosion cycle, which is also extended by the effect of roots through delaying entrainment rates. In contrast, more erodible soils fail at earlier stages even if root-reinforced, which reduces the time scale of the entrainment phase and thus the period along which roots reduce shear stresses onto the soil. As a consequence, this positive feedback enlarge further the difference of erosion rates between sandy and loam soils in presence of woody roots.

Floodplain heterogeneity on bank erosion and river migration
The angle between past and present river channels, the stratification of the floodplain, and the regulation of the water levels, created the conditions for wide-ranging erosion rates along relatively short distances. The resulting scale of the embayments is not common in rivers where channels normally follow abandoned paths Constantine et al., 2010) or cut previous paths at large angles (Hooke, 1995;Slingerland and Smith, 2004). Yet, as expected, large embayments also arose in other reaches of the Meuse River, for instance, Previous studies proposed flow patterns and their interaction with bank shape as the main factors for the formation of embayments (Hackney et al., 2015). However, the appearance and growth of the embayments in the Meuse River are dominated by floodplain stratification, under relatively uniform and highly delimited loads induced by ship waves. Therefore, in cases where flow velocities are not affected by stage regulation, bank stratigraphy and composition may also have a significant role in shaping banklines into large embayments. This hypothesis is supported by the fact that the bank bathymetric surveys performed by Hackney et al. (2015) in the Mekong River presented clear strata, and that critical shear stresses ranged significantly from the weakest to the strongest soil samples, with respect to acting shear stresses.
This work indicates that structured along-channel variability in lithology defines erosion patterns within a river reach, at a scale which has not been identified before. Previous works have shown the control that floodplains with heterogeneous erosional resistance exert on meander migration, through both field observations (Hudson and Kesel, 2000) and numerical models . Further numerical analyses and field cases showed that meander planform complexity increases with floodplain soil heterogeneity Vermeulen et al., 2014). Konsoer et al. (2016b) found that vertical heterogeneity in bank resistance significantly changes bank erosion rates and mechanisms in the Wabash River, in agreement with Motta et al.'s (2014) long-term meander simulations, which included a physically-based erosion model for stratified banks.
This work highlights the importance of floodplain formation and related channel dynamics for later morphological developments. The complexity of the processes involved during floodplain formation (Nanson and Croke, 1992;Kleinhans, 2010;Kleinhans et al., 2018) limits the use of uniform bank erosion coefficients for predictive purposes, even considering complex hydrodynamics . Moreover, Schwendel et al. (2015) showed significant changes in migration rates and sinuosity of a large meandering sand-bed river driven by clay bodies, presumably caused by pedogenic processes.  demonstrated that numerical models that account for floodplain formation, including key geomorphic units as scroll bars and oxbow lakes, achieve realistic meander planforms thanks to heterogeneous bank resistances. The morphological features within a river reach, even when developing in the short term (<10 years), may have an impact for further developments and latter affect landscape formation, stressing the relevance of process-based numerical modelling.

Conclusions
Distinct bankline patterns, presenting several oblique embayments with different angulations, arose along two reaches of the Meuse River after removing bank protections. We analysed these patterns by considering the processes and factors affecting bank erosion rates in the study area. The results show that floodplain heterogeneity controls the bankline irregularities in magnitude and orientation, in combination with ship waves hitting the bank at regulated water levels. Past river alignments and floodplain stratigraphy indicate that scroll bar depositions modulate the erodibility along the channel. Different strata thicknesses and orientations with respect to the river channel can explain the different embayment lengths and the stretches with parallel bankline retreat.
On the bank strength side, the effectiveness of isolated tree roots on reducing bank erosion rates depends on soil characteristics, which define a primary control over erosion rates. Locations with sandy deposits vanished the contribution of mature poplar trees over yearly intervals, whereas loam layers seem to allow tree roots to significantly increase the resistance against entrainment and reduce bank erosion rates. On the load side, floods induce relatively low bed shear stresses during short periods, whereas ship waves regularly exert similar to higher loads. The difference in frequency and duration between these loads, and the fact that water levels are strongly regulated, create the conditions for ship waves to form a long mild terrace across the banks.
This research highlights the importance of underlying mechanisms that act at large spatial scales and control future responses of eroding riverbanks. This supports the idea that process-based long-term meandering models should account for processes that form floodplains to improve the representation of natural planforms and their eventual predictive accuracy, advancing statistical approaches. Furthermore, shortterm and process-based bank erosion models applied to navigable rivers need to consider heterogeneous floodplain properties, but also have the challenge to represent the magnitude of loads exerted by ship waves over different phases of bank evolution.