Ripples formed and migrated downstream in the gaps between patches. For example, point A in Figure 2A (red symbols) was located between patches. At this point, the bed elevation fluctuated irregularly over Δ z = ± 0.4 cm for the first 3 h, then fluctuated periodically with an amplitude of 0.7 cm and a period of 80 min. In contrast, directly behind a patch (point B, blue symbols in Figure 2A ), sediment was continuously deposited over time, resulting in a steady increase in bed elevation. This observation aligned with previous studies that also reported sediment deposition downstream of a vegetation patch (e.g., Ortiz et al., 2013; Hu et al., 2018).

The bedload transport rate was highest in the bare channel and decreased with increasing ϕ p ( Figure 2B ). Given the limited number of cases, a simple linear decrease would be a reasonable inference. However, when considering the 2D (red) and 3D (yellow) cases in isolation, the decrease in q s was steeper for ϕ p > 0.3, suggesting a threshold behavior, which we explore further in the context of the velocity and turbulence trends. The trend of decreasing bedload transport with increasing vegetation coverage contrasts with emergent patchy canopies, where bedload transport rate increases as the bed area occupied by vegetation ( ϕ p ) increases (Shan et al., 2020). In emergent canopies, flow is laterally deflected around the patches, causing velocity and turbulence to increase in the bare areas between vegetation, which in turn enhances the bedload transport. Conversely, in submerged canopies, as the bed area occupied by vegetation increases, more flow passes over the top of individual patches rather than around them, reducing the near-bed velocity ( Figure 3C ) and consequently diminishing q s . Finally, for the same ϕ p , shifting from a distribution of finite patches (3D pattern, Figure 1C ) to channel-spanning patches (2D pattern, Figure 1C ) resulted in a lower bedload transport rate (compare red and orange dots at same ϕ p in Figure 2B ). The variation in q s was influenced by both mean and turbulent flow characteristics, which will be discussed in the next section.

Near-bed flow characteristics differed between spatial configurations. First, consider case 3 with a patch area density of ϕ p = 0.3. With large lateral gaps between individual patches, near-bed flow primarily moved through the gaps, with significantly lower velocity ( u ¯ / U 0 ≈ 0) and lower TKE ( k t / U 0 2 ≈ 0.02) downstream of individual patches (green shading in Case 3 in Figure 3A (2), Figure 3B (2); Supplementary Figure S1 ) compared to the regions between the patches (brown shading in Figure 3A (2), Figure 3B (2); Supplementary Figure S1 ). The regions of low hydrodynamic intensity behind the patches facilitated sediment deposition (point B in Figure 2A ). In the regions in between patches, higher mean velocity and TKE promoted ripple formation and sediment migration (point A in Figure 2A ). As the lateral spacing between patches decreased, the velocity difference between the gaps and wakes also decreased (compare cases 3 and 5 in Figure 3A ).

The spatially averaged near-bed velocity, U n b , was highest in the bare channel (blue marker, Figure 3C ) and decreased with increasing patch area density ϕ p . Patch configuration also influenced velocity at the same ϕ p : the 2D channel-spanning configurations (red markers in Figure 3C ) consistently exhibited lower mean velocity than 3D configurations (orange markers in Figure 3C ) due to total blockage when patches spanned the channel width as a continuous group (2D pattern in Figure 1C ). For ϕ p = 0.56 in a channel-spanning configuration, most of the flow skimmed over the canopy, reducing near-bed velocity to nearly zero ( Figure 3A , Case 9). In this case, the ratio of surface velocity ( 〈 u ¯ 〉 m a x ) to near-bed velocity ( 〈 u ¯ 〉 m i n ) was 0.048, which met the criteria for skimming flow, 〈 u ¯ 〉 m i n / 〈 u ¯ 〉 m a x < 0.05, as defined in Folkard (2011). ϕ p = 0.56 also met the criteria ϕ p > 0.4 for skimming flow (Wolfe and Nickling, 1993).

In contrast to near-bed velocity, which monotonically decreased with increasing ϕ p ( Figure 3C ), near-bed turbulence exhibited a nonlinear trend with ϕ p , with a peak near ϕ p = 0.3 (circles in Figure 3D ). Given the limited number of cases, it was possible that the peak was broad and spanned ϕ p = 0.3 and 0.4. Specifically, considering all cases at the same ϕ p , the average TKE for ϕ p = 0.3 and 0.4 were 0.045 ± 0.028 and 0.033 ± 0.005, respectively ( ± 95% confidence level in a t-test), and the probability that the confidence intervals of these two groups overlap was approximately 32%. Figure 3D is supplemented with measurements from Park and Nepf (2025, triangles in Figure 3D ), who conducted experiments with the same vegetation models but over a fixed bed. For similar canopy conditions, a lower TKE was observed over the mobile bed (circles), which likely reflected the loss of turbulent energy transferred to sediment movement. Despite the difference associated with bed mobility, both studies indicate a nonlinear influence of patch area density and near-bed TKE, with a maximum near ϕ p = 0.3 ( Figure 3D ).

For a bare bed (blue circle in Figure 3D ), near-bed turbulence was solely derived from bed-shear production ( P s ). The introduction of patchy vegetation both changed the near-bed shear production and created new sources of turbulence, resulting in higher near-bed TKE (red and orange markers in Figure 3D ) compared to the bare bed (blue circle). Specifically, wake production ( P w ) occurred in the canopy, and vertical turbulent transport ( T t ) and dispersive transport ( T d ) transferred turbulence generated in the shear-layer at the canopy top to the near-bed region (see details in Supplementary Data Sheet 1 - Section 2). A notable exception was the canopies classified as skimming flow (green markers in Figure 3D ), for which near-bed turbulence was lower than the bare bed. The nonlinear response of near-bed TKE to patch area density can be understood through changes in the turbulence budget.

When ϕ p was low ( ϕ p = 0.13), near-bed velocity was high ( Figure 3C ), resulting in high wake production ( P w ), which contributed the most (38%) to the near-bed turbulence budget (red bar in Figure 3(G)-1 ). As patch density increased to ϕ p = 0.3, near-bed velocity diminished and became more laterally uniform, reducing the contributions from both bed-shear production and wake production (orange and red bars in Figure 3G (2)). At this density, turbulence produced at the top of the canopy and carried toward the bed by turbulent transport ( T t ) made the greatest contribution to near-bed TKE ( Figures 3(G)-2 and 3(G)-3 ).

For the patch area density producing the peak TKE ( ϕ p = 0.30), the magnitude of near-bed turbulence depended on spatial configuration, with lower 〈 k t 〉 / U 0 2 = 0.03 for distributed 3D patches (orange markers in Figure 3D ) compared to 〈 k t 〉 / U 0 2 = 0.06 for channel-spanning 2D patches (red marker in Figure 3D ). The higher turbulence in the channel-spanning configuration was attributed to strong shear in the canopy shear layer, which contributed to the near-bed turbulence through a higher magnitude of vertical turbulent transport ( T t , blue bars in Figure 3G ).

As patch area density increased to the point of skimming flow ( ϕ p = 0.56), the shear at the top (green squares in Figure 3E ) and thus the canopy-shear production reached its highest value. However, the turbulence produced at the canopy top did not penetrate to the bed, as indicated by the near-zero T t below z / h = 0.4 (green square in Figure 3F ). The high canopy drag inhibited turbulent transport toward the bed. In addition, under skimming flow conditions, the near-bed velocity was both small and uniform in both lateral and vertical directions ( Figures 3(A)-4 , 3E ), resulting in minimal bed-shear and wake production. Consequently, the minimum near-bed TKE was observed among all cases (green marker in Figure 3D ). Finally, dispersive transport did not make a significant contribution to near-bed turbulence in any of the cases, consistently remaining below 15% (green bars in Figure 3G ).

The variation in velocity and TKE with ϕ p explained the variation of q s with ϕ p . Recall from Figure 2B that q s decreased as ϕ p increased, but the trend was milder for ϕ p < 0.3 and stronger for ϕ p > 0.3. For ϕ p < 0.3, as ϕ p increased, near-bed velocity decreased, but TKE increased ( Figures 3C, D ). These trends had opposing impacts on bedload transport, resulting in a mild decrease in q s with increasing ϕ p ( Figure 2B ). However, for ϕ p > 0.3, both velocity and TKE decreased with increasing ϕ p ( Figures 3C, D ), causing a sharp decrease in q s ( Figure 2B ). Under skimming flow conditions at the highest ϕ p , both mean velocity and TKE reached their lowest value across all experimental cases, resulting in the smallest q s . This highlights the importance of considering both time-mean velocity and turbulence when interpreting bedload transport within submerged canopies.

The drop in bedload transport near ϕ p = 0.3 can be linked to a positive feedback mechanism known to influence the vegetated landscapes. Specifically, in seagrass meadows, which predominantly experience submerged conditions, a patch density around ϕ p = 0.4 is commonly observed in fragmented meadows (e.g., Luhar et al., 2008; Schaefer et al., 2025). According to percolation theory, when ϕ p < 0.4, bare regions can form channels associated with higher velocity, promoting sediment mobility and inhibiting vegetation re-establishment. This creates a positive feedback loop that sustains the channels and the patchy meadow landscape. For example, Schaefer et al. (2025) used satellite images to show that a patchy seagrass meadow persisted for over 50 years, with individual patches shifting position while maintaining ϕ p ≈ 0.4. In contrast, when ϕ p > 0.4, continuous channels are less likely to form, and velocity in bare regions is not elevated, eliminating the positive feedback. Consequently, when ϕ p > 0.4, patches can expand into the bare regions and eventually form a continuous meadow ( ϕ p = 1).

Our results revealed a similar threshold between ϕ p = 0.3 and 0.4. Specifically, when ϕ p ≤ 0.3, the bare regions between patches exhibited higher near-bed velocity and greater bedload transport, potentially inhibiting vegetation regrowth and maintaining the existing vegetation distribution. However, when ϕ p ≥ 0.4, both mean velocity and TKE sharply decreased, causing bedload transport to approach zero as ϕ p increased ( Figure 2B ). The reduction in bedload transport would stabilize the channel bed, protecting the existing vegetation distribution and potentially facilitating vegetation expansion into the bare regions. Thus, our study aligns with previous studies suggesting that a threshold near ϕ p = 0.3 to 0.4 plays a key role in the evolution of vegetated landscapes.

The measured bedload transport was compared with the models introduced in Section 2, which were based on bed shear stress (Equation 1), near-bed TKE (Equation 4), and a combination of both (Equation 5). The bed drag coefficient C f , n b ( = 0.005 ± 0.001) and ω ( = 0.11) were computed from measurements over a bare bed, as described in Section 2. For the bare bed condition, both the τ and k models provided good prediction (blue triangles in Figures 4A, B ), but the τ model significantly underestimated q s in the vegetated channels, with errors reaching up to two orders of magnitude ( Figure 4A ). This was because the time-averaged bed shear stress ( τ ) did not account for vegetation-generated turbulence, which can play a critical role in bedload transport within vegetated channels (Tinoco and Coco, 2016, 2018; Zhao and Nepf, 2021). The k model provided a better prediction for emergent canopies, compared to the τ model (see red markers in Figures 4A, B ), but it performed poorly for submerged canopies (blue circles and triangle markers in Figure 4B ). For submerged canopies, turbulence transported vertically from the top of the canopy contributed significantly to near-bed turbulence (Zhao and Nepf, 2024), leading to a weak correlation between near-bed velocity and TKE, with a correlation coefficient of just 0.03 ( Supplementary Figure S2 ). In this case, neither the bed shear stress estimated from time-averaged velocity nor the TKE alone could adequately describe the flow characteristics influencing bedload transport.

The hybrid model defined in Equations 5 and 6 provided a closer fit to observations for both bare and vegetated channels and, importantly, for both emergent and submerged canopies ( Figure 4C ). However, in this study, using patchy submerged vegetation models, the value of β in Equation 5 was found to be β = 7, which was lower than the previously reported β = 10 in studies using both emergent and submerged cylinder arrays. The lower β can be explained by the difference in canopy structure. Specifically, in patchy canopies, stronger and more pronounced regions of horizontal shear formed around individual patches compared to individual cylinders, which elevated turbulence levels and resulted in a higher k n b / U n b 2 compared to the cylinder arrays in previous studies ( Figure 4D ). The lower β in the patchy canopy also implied that, under the same near-bed velocity and TKE conditions, the patchy canopy produced lower bedload transport compared to cylinder arrays. This could be attributed to the larger bedform development within the patchy canopy. While cylinder arrays generate localized scour holes around individual cylinders, the patchy canopy formed relatively large deposition regions downstream of individual patches ( Figure 2A ). These deposition zones did not contribute to channel-scale bedload transport, leading to a lower q s , * compared to cylinder arrays.

The present study considered a single sediment size ( d s = 0.35 mm), and the results would likely be influenced by sediment characteristics that affect their response to turbulence. For example, the response timescale for a sand grain with d s = 0.35 mm is t p = ρ s d s 2 / 18 μ = 0.018 s, where μ is the dynamic viscosity. The near-bed turbulence timescale ( t f ), defined as the integral time scale, varied from 0.18 to 0.26 s across the experimental cases. These time scales resulted in Stokes numbers ( S t = t p / t f ) below 0.1, indicating that the sediment grains could respond to the turbulent fluctuations (Brandon and Aggarwal, 2001). However, for d s > 1.3 mm, the sediment response time would be longer than near-bed turbulence time scales ( S t > 1), meaning individual grains would not respond to turbulent fluctuations, and sediment mobility would be governed solely by the time-mean bed stress. A similar time-scale control has been observed in biocohesive sediment. Specifically, biofilms can bind individual sediment grains together, increasing their response time (Malarkey et al., 2015). As observed by Deitrick et al. (2023), when the time scale of the vegetation-generated turbulence is shorter than the time needed to break biocohesion between sediment grains, the turbulence has no impact on sediment transport, and transport is best described by the mean shear stress.