Statistical Analysis of Bed Load Transport over an Armored Bed Layer with Cluster Microforms

: River engineers have long been challenged by the need to predict sediment transport, especially over armored riverbeds. This study investigates the statistical properties of bed load transport over an armored bed layer with cluster microforms in laboratory experiments. Particle clusters on the sediment bed were formed by widely graded particles under constant ﬂow. A series of key kinematic parameters computed from particle trajectories recorded by a digital camcorder, including mean squared particle displacement (MSD), particle number activity, particle velocities, step length, and rest period, were analyzed. The scaling growth of the MSD with time showed that the particle di ﬀ usion regime was superdi ﬀ usive at small time scales, but became subdi ﬀ usive at larger time scales. The particle number activity follows a negative binomial distribution, and the probability distributions of streamwise and transverse particle velocities displayed heavy asymptotic tails, which indicates the particle clusters might exert a dual impact on bed load transport: some particles are accelerated in the preferential paths between particle clusters, while others were obstructed by the particle clusters. In addition, the bed load di ﬀ usion regime varied with observation time scales. The ﬁndings of this study can gain insight into the bed load transport processes over armored riverbeds. time scales associated with the transition of particle di ﬀ usion regimes and those with the changes in the statistical characteristics of bed load particle motions. Further investigation is needed in unraveling their connections.


Introduction
In nature, a gravel stream bed that has a wide particle size distribution often develops an armored surface layer with various cluster microforms [1][2][3]. The term "cluster" describes a discrete, self-organized grouping of sediment particles that sit atop the surrounding bed surface [4]. Over the last two decades, river engineers have been challenged by the need to predict sandy sediment transport over armored gravel beds. Such phenomena are often observed in streams after dam removal or wildfires [5]. Despite the substantial efforts, including hiding-exposing correction [6] and surface-based transport mechanisms [7], our ability to estimate bed load flux on mixed-size gravel beds remains unsatisfactory [8]. A major obstacle is to quantify the impact of particle clusters on bed load transport.
Bed load transport is a highly stochastic process, which involves an enormous number of particle motions that can be viewed as a series of probabilistic events. It remains a challenge to explore the particle motions at the grain scale even in a well-controlled flume experiment due to complex flow-particle, particle-particle and particle-bed interactions. Fortunately, recent advances in image-processing technology allow us to track the trajectories of multiple sediment particles simultaneously, which provides a means to capture the details of instantaneous bed load particle motions [9][10][11][12] rather than the average quantities only.
process. The results from our experiment may provide data to calibrate the key parameters of these models and gain insight into the statistical properties of bed load transport over highly nonuniform sediment beds.

Theoretical Backgrounds
Nikora et al. [39] established a theoretical framework, which generalized the bed load transport as a particle diffusion process. Bed load transport can be viewed as an ensemble of sediment particles, each of which moves episodically through discrete steps separated by resting periods. Although the forms of individual particle motions may differ (i.e., rolling, sliding, and saltating), the bed load transport can be conceptualized as a random walk with two distinct states, namely the moving state and the resting state. In other words, we do not discriminate between different forms of particle motions and simply classify all active particles as "moving" particles. The diffusion regime of bed load transport can be determined by the temporal variation in displacements of bed load particles. The MSD, denoted by σ 2 , may be adopted to define the diffusivity of the bed load transport process, which can be expressed as [36]: where x(t) is the displacement of an individual particle at time t, and the angle bracket denotes an ensemble average over all particles. The diffusion regime of bed load transport can be identified as normal (Fickian) or anomalous (non-Fickian) according to the scaling growth of the MSD of particle displacements with time t [15,39]: where γ is a power-law exponent of t. For normal particle diffusion, the exponent γ = 1 and σ 2 grows linearly with time; for subdiffusion, 0 < γ < 1; for superdiffusion, γ > 1; and for ballistic transport, γ = 2 and σ 2 grows the same as particles moving in different directions with no pauses or direction changes [36]. Theoretically speaking, the diffusion regime can be fully described by the collective behavior of bed load particles at the grain scale. Although some conceptual models, e.g., the symmetric/asymmetric random walks [36] and Langevin-Fokker Planck equations [15], have been proposed to quantify the connection between the particle diffusion regime and individual particle motions, the details have not been fully revealed. Normal diffusion is a corollary of the Central Limit Theorem (CLT), which states that the sum of a large number of statistically independent events will be Gaussian distributed, provided that the distribution of the individual events is not too broad. In case of bed load transport, the CLT requires that both the step length and rest period have converging first and second moments for normal diffusion to occur. For example, if both the step length and the rest period have the exponential probability distributions, the transport system is normally diffusive [36]. The exponential distributions for the step length s and the rest period t can be expressed in terms of the mean step length S and the mean rest period T as: In normal diffusion, the probability in which a particle is found at a specific location decreases exponentially with the distance away from its original location. If the distributions of the step length and/or the rest period have thicker (heavier) tails than do the exponential distributions, anomalous diffusion may result. Weeks et al. [36] reported strong dependence of the diffusion regimes on the heavy-tailed distributions of the step length and the rest period. In their model, both the step length and the rest period are described by power-law distributions, which can be expressed as: where µ, ν are the exponents characterizing the power-law decay rates. When µ, ν < 3, the distributions of the step length and the rest period are broad, implying an unusually large number of long particle steps and long particle rests than predicted by Gaussian models. Other heavy-tailed PDFs, such as the Pareto distribution [15] and the Mittag-Leffle (M-L) distribution [40], were also shown to significantly affect the diffusion regimes of bed load transport [11,15]. Particle number activity and particle velocity are also important stochastic parameters, the PDFs of which differ between normal and anomalous diffusion. The two components of bed load particle motions are key elements in describing the sediment flux and the dispersal of particles during transport. Instantaneous or time-averaged gradients in the particle activity affect particle flux, and the second moment of the particle velocity distribution influences the rate of diffusion of bed load particles in the presence of a particle activity gradient [38]. Roseberry et al. [13] defined the particle number activity as the number of moving particles per unit area per unit time. In the present study, the particle number activity is represented by the number of moving particles within the observation window between two consecutive video snapshots, which is essentially equivalent to that defined by Roseberry et al. [13] because the window size and the frame rate remain unchanged during our experiment. The particle number activity is a key element in describing the diffusive part of the bed load flux [38]. Ancey et al. [22] derived an analytical solution for particle number activity under equilibrium transport, which takes the form of negative binomial distribution as: where n is particle number activity, r, p are parameters related to particle entrainment and deposition probabilities, and Γ() denotes a gamma distribution. Conventionally, the streamwise particle velocity was modeled by the exponential or Gaussian distributions [15,41]. For example, the exponential distribution for streamwise particle velocity can be expressed as [10,13,17,41]: where V x is the mean particle velocity at equilibrium transport. Some researchers proposed that the PDFs of the streamwise particle velocity may be heavy-tailed in anomalous diffusion. For example, Bradley et al. [42] modeled the streamwise particle velocity in anomalous diffusion system with an α-stable distribution with a power-law asymptotic tail. The characteristic exponent α indicates the deviation from normal diffusion: the α-stable distribution is Gaussian if α = 2, and is heavy-tailed otherwise. Readers may refer to Bradley et al. [42] for more description on the α-stable distributions. The cross-stream particle velocity was proposed to be modeled by the Gaussian distribution, which can be expressed as [10,17]: where V y is the cross-stream particle velocity, and δV y characterizes the order of magnitude of the cross-stream velocity fluctuation. It is worth noting that the probability distributions reviewed in this section were proposed for bed load transport over uniform or well-sorted (narrowly graded) sediment beds. For bed load transport over nonuniform sediment beds, it is of importance to find a proper distribution that may well describe the key parameters of bed load particle motions and analyze their connection with the diffusion regimes of the bed load transport.

Experimental Setup
The experiment was conducted in a rectangular flume at China Institute of Water Resources and Hydropower Research (IWHR). The IWHR flume was 50 m long, 1 m wide, and 1.2 m deep with fixed glass walls on both sides of the flume, and the bed slope was 0.004 ( Figure 1). A fixed-bed section (27 m in length) was placed after the flume entrance, which was packed with pebbles and cobbles in order to stabilize flow and dissipate flow energy. The test section (20 m in length) was placed between the fixed-bed section and the flume outlet, which was paved with widely graded sediment ranging from sand (0.25 mm) to pebbles (20 mm) ( Table 1). The particle size distribution of the sediment bed was selected according to the sediment data at the Lujiahe Site in the middle Yangtze River, China, where strong bed armoring processes were observed. A bed load sediment trap was installed at the flume outlet to collect bed load sediment samples. A digital camcorder was mounted on top of the test section to record the motions of bed load particles at 10 frames per second with a resolution of 1280 × 720 pixels. The monitoring area was a 366 mm (streamwise) × 206 mm (cross-stream) region of the bed surface near the downstream end of the test section. A sheet of glass was placed on the water surface to reduce image distortion due to water surface oscillation.

Experimental Setup
The experiment was conducted in a rectangular flume at China Institute of Water Resources and Hydropower Research (IWHR). The IWHR flume was 50 m long, 1 m wide, and 1.2 m deep with fixed glass walls on both sides of the flume, and the bed slope was 0.004 ( Figure 1). A fixed-bed section (27 m in length) was placed after the flume entrance, which was packed with pebbles and cobbles in order to stabilize flow and dissipate flow energy. The test section (20 m in length) was placed between the fixed-bed section and the flume outlet, which was paved with widely graded sediment ranging from sand (0.25 mm) to pebbles (20 mm) ( Table 1). The particle size distribution of the sediment bed was selected according to the sediment data at the Lujiahe Site in the middle Yangtze River, China, where strong bed armoring processes were observed. A bed load sediment trap was installed at the flume outlet to collect bed load sediment samples. A digital camcorder was mounted on top of the test section to record the motions of bed load particles at 10 frames per second with a resolution of 1280 × 720 pixels. The monitoring area was a 366 mm (streamwise) × 206 mm (cross-stream) region of the bed surface near the downstream end of the test section. A sheet of glass was placed on the water surface to reduce image distortion due to water surface oscillation.   Before the experiment started, the sediment bed was leveled and worked by uniform water flow to form an armored sediment bed with cluster microforms. The discharge was initially maintained at 80 L/s and then increased in a stepwise manner to 170 L/s. At each flow rate, the flow was maintained constant for sufficiently long time until the flume bed reached an equilibrium state. Moving grains were collected by the sediment trap installed at the flume outlet and continuously fed at the flume inlet. A quasi-steady armored bed layer was formed over the mobile section of the flume after 78 h, as the bed load transport rate approached a small and stable flux rate between 0.74 g/s and 1.05 g/s with an average value of 0.84 g/s for 9.3 h.
After the armored bed layer was formed, the tracer experiment was started with a constant flow rate of 130 L/s. The hydraulic parameters of the experiment were summarized in Table 2. The 3 mm black grains were used as tracer particles and fed at the flume entrance. The motions of the tracer   Before the experiment started, the sediment bed was leveled and worked by uniform water flow to form an armored sediment bed with cluster microforms. The discharge was initially maintained at 80 L/s and then increased in a stepwise manner to 170 L/s. At each flow rate, the flow was maintained constant for sufficiently long time until the flume bed reached an equilibrium state. Moving grains were collected by the sediment trap installed at the flume outlet and continuously fed at the flume inlet. A quasi-steady armored bed layer was formed over the mobile section of the flume after 78 h, as the bed load transport rate approached a small and stable flux rate between 0.74 g/s and 1.05 g/s with an average value of 0.84 g/s for 9.3 h.
After the armored bed layer was formed, the tracer experiment was started with a constant flow rate of 130 L/s. The hydraulic parameters of the experiment were summarized in Table 2. The 3 mm black grains were used as tracer particles and fed at the flume entrance. The motions of the tracer particles were filmed by the 10 Hz digital camcorder when they passed the monitoring area, which can be processed using image processing techniques introduced in the following section.

Image Processing
The trajectory of a tracer particle can be obtained from a sequence of image frames extracted from video clips. There are two types of popular automated techniques in tracking tracer particles, namely the "shape-matching" approach and the "pattern-matching" approach [43,44]. However, both approaches may encounter difficulty if the tracer particles experience unusually long particle steps [29,39,42]. In our study, the centroids of the tracer particles in each image frame were manually marked with ImageJ program to avoid computer misjudgment. ImageJ is a public domain, Java-based image-processing software developed by the United States National Institutes of Health that has been widely used to track particle trajectories in sediment transport studies [11,13]. The particle trajectory can be automatically generated from the marked locations of tracer particles and converted into a series of streamwise and transverse coordinates, denoted by (x i , y i ), that are functions of time t i . The streamwise and transverse particle velocities, V x and V y , can be readily computed as: Considering small image distortion may exist, particle steps that are less than one pixel were excluded from the dataset. Because it is difficult to distinguish between the rolling, sliding and saltating modes of bed load particles on plan-view images, we define the step length as the distance a particle covers between two consecutive "resting" states. A particle was considered to be "resting" if the displacements in both the x-and y-directions in two sequential image frames were less than one pixel; otherwise, it was considered to be "moving". Similarly, the rest period is defined as the period between two consecutive particle steps.

Experimental Results
In this study, we tracked the 3 mm tracer particles passing through the monitoring area in a 5 min video clip, in which 1267 particle steps containing over 300,000 centroid locations can be identified. We analyzed the temporal variation in the MSD and the probability distributions of the particle number activity, the particle velocity, the step length, and the rest period in Sections 4.1-4.4, respectively. Fan et al. [45] argued that all tracers must be tracked because the loss of tracked particles may result in apparent but deceptive subdiffusion. In our experimental dataset, only 1% of the tracers were not tracked during the study period.

MSD Growth and Particle Diffusion Regimes
The variation in the MSD of the tracer particles with time was plotted in Figure 2. When fitted by the power-law function in the form of Equation (1), the relationship of the MSD versus time showed different power-law exponents at different time scales: the left part of the fitted curve has a power-law exponent α > 1, while the central and the right parts of curve have power-law exponents α < 1. This scale-dependent behavior is consistent with the conceptual model with time scales proposed by Nikora et al. [39], who proposed that bed load particle motions comprised of three ranges of temporal scales, namely the local range, the intermediate range, and the global range, which were associated with different particle diffusion regimes. In the present study, the particle diffusion regime was superdiffusive in the range of t < 50 s, and was subdiffusive at larger time scales. Fan et al. [15] suggested that the anomalous diffusion regime in bed load transport is resulted from the heavy-tailed distributions of the kinematic variables of particle motions. In addition, the time scales at which the transitions between particle diffusion regimes occur are apparently larger than that proposed by Nikora et al. [39] and Fan et al. [15], which might also be due to the thick tails in the probability distribution of the rest period. Therefore, the key kinematic parameters, including particle number activity, particle velocities, step length, and rest period were analyzed in the following sections. ranges of temporal scales, namely the local range, the intermediate range, and the global range, which were associated with different particle diffusion regimes. In the present study, the particle diffusion regime was superdiffusive in the range of t < 50 s, and was subdiffusive at larger time scales. Fan et al. [15] suggested that the anomalous diffusion regime in bed load transport is resulted from the heavy-tailed distributions of the kinematic variables of particle motions. In addition, the time scales at which the transitions between particle diffusion regimes occur are apparently larger than that proposed by Nikora et al. [39] and Fan et al. [15], which might also be due to the thick tails in the probability distribution of the rest period. Therefore, the key kinematic parameters, including particle number activity, particle velocities, step length, and rest period were analyzed in the following sections.

Particle Number Activity
The probability distribution of the particle number activity was plotted in Figure 3. Ancey et al. [22] considered bed load transport as a two-state Markov process and derived that the particle number activity followed the negative binomial distribution if collective entrainment (entrainment of stationary grains by active ones) is considered or the Poisson distribution otherwise. Therefore, the Poisson distribution and the negative binomial distribution were plotted in the figure for comparison. Figure 3 shows that the measured data significantly deviate from the Poisson distribution, but is reasonably well described by the negative binomial distribution, which indicates that the particle interaction could potentially affect the probability distribution of the particle number activity. However, particle collisions between bed load particles were not frequently seen in our experiment due to relatively low bed load transport rates. Roseberry et al. [13] suggested that the discrepancy between the distribution of the actual particle number activity and the Poisson distribution might also be attributed to the turbulence structure induced by bed configuration, leading to preferential selection or exclusion of particle entrainment and deposition. In our experiment, cluster formation and breakdown were occasionally observed in the armored bed layer; consequently, some tracer particles were trapped inside the particle clusters, whereas others were preferentially transported via the paths between the particle clusters. Therefore, the preferential particle transport and trapping by the particle clusters are likely to be responsible for the discrepancy between the two probability distributions.

Particle Number Activity
The probability distribution of the particle number activity was plotted in Figure 3. Ancey et al. [22] considered bed load transport as a two-state Markov process and derived that the particle number activity followed the negative binomial distribution if collective entrainment (entrainment of stationary grains by active ones) is considered or the Poisson distribution otherwise. Therefore, the Poisson distribution and the negative binomial distribution were plotted in the figure for comparison. Figure 3 shows that the measured data significantly deviate from the Poisson distribution, but is reasonably well described by the negative binomial distribution, which indicates that the particle interaction could potentially affect the probability distribution of the particle number activity. However, particle collisions between bed load particles were not frequently seen in our experiment due to relatively low bed load transport rates. Roseberry et al. [13] suggested that the discrepancy between the distribution of the actual particle number activity and the Poisson distribution might also be attributed to the turbulence structure induced by bed configuration, leading to preferential selection or exclusion of particle entrainment and deposition. In our experiment, cluster formation and breakdown were occasionally observed in the armored bed layer; consequently, some tracer particles were trapped inside the particle clusters, whereas others were preferentially transported via the paths between the particle clusters. Therefore, the preferential particle transport and trapping by the particle clusters are likely to be responsible for the discrepancy between the two probability distributions.

Particle Velocity
The probability distribution of the streamwise particle velocity, Vx, was plotted in Figure 4. The probability distributions of the particle velocity were investigated in several previous experiments. Lajeunesse et al. [10], Roseberry et al. [13] and Shim and Duan [17]neglected negative streamwise velocities in their experiment data and found that the rest of the data were well fitted by the exponential distribution, whereas Martin et al. [11], Ancey and Heyman [41]and Heyman et al. [12] suggested that the particle velocity followed the Gaussian distribution. Fan et al. [20] used a hybrid mechanistic-stochastic model to analyze the bed load transport under low transport rates and found that the positive and negative streamwise velocities decayed exponentially at different rates. In our experiment, negative velocities accounted for almost 10% of the total Vx measurements and therefore cannot be neglected. Therefore, we adopted the (thin-tailed) Gaussian distribution with exponential tails and the (heavy-tailed) α-stable distribution with power-law tails to describe the two-sided Vx distribution. Figure 4 shows that the probability distribution of the streamwise particle velocity significantly deviates from the Gaussian distribution as Vx→+∞ and is well described by the α-stable distribution. The thick tails in the α-stable distribution imply that there is a higher fraction of highvelocity bed load particles over the armored bed layer than predicted by the thin-tailed distributions. The left tail in the α-stable distribution on the left indicates a noticeable amount of grains moved in the upstream direction, which has not been reported in other experiments.

Particle Velocity
The probability distribution of the streamwise particle velocity, Vx, was plotted in Figure 4. The probability distributions of the particle velocity were investigated in several previous experiments. Lajeunesse et al. [10], Roseberry et al. [13] and Shim and Duan [17] neglected negative streamwise velocities in their experiment data and found that the rest of the data were well fitted by the exponential distribution, whereas Martin et al. [11], Ancey and Heyman [41] and Heyman et al. [12] suggested that the particle velocity followed the Gaussian distribution. Fan et al. [20] used a hybrid mechanistic-stochastic model to analyze the bed load transport under low transport rates and found that the positive and negative streamwise velocities decayed exponentially at different rates. In our experiment, negative velocities accounted for almost 10% of the total V x measurements and therefore cannot be neglected. Therefore, we adopted the (thin-tailed) Gaussian distribution with exponential tails and the (heavy-tailed) α-stable distribution with power-law tails to describe the two-sided V x distribution. Figure 4 shows that the probability distribution of the streamwise particle velocity significantly deviates from the Gaussian distribution as V x →+∞ and is well described by the α-stable distribution. The thick tails in the α-stable distribution imply that there is a higher fraction of high-velocity bed load particles over the armored bed layer than predicted by the thin-tailed distributions. The left tail in the α-stable distribution on the left indicates a noticeable amount of grains moved in the upstream direction, which has not been reported in other experiments. tails and the (heavy-tailed) α-stable distribution with power-law tails to describe the two-sided Vx distribution. Figure 4 shows that the probability distribution of the streamwise particle velocity significantly deviates from the Gaussian distribution as Vx→+∞ and is well described by the α-stable distribution. The thick tails in the α-stable distribution imply that there is a higher fraction of highvelocity bed load particles over the armored bed layer than predicted by the thin-tailed distributions. The left tail in the α-stable distribution on the left indicates a noticeable amount of grains moved in the upstream direction, which has not been reported in other experiments.  The probability distribution of the transverse particle velocity, V y , was plotted in Figure 5. Similarly, V y shows a symmetric bell-shaped distribution that has a narrow central peak and heavy tails on both sides that are much thicker than those of the Gaussian distribution. Thus, we proposed that the α-stable distribution also provides a better fit for V y with proper parameters. Our experimental results contradict the widely accepted view that the transverse particle velocity follows the thin-tailed exponential or Gaussian distributions based on previous experiments [10,13,17]. After comparing with the previous experiments, we attribute this inconsistency to the armored gravel bed configuration with cluster microforms, which was not commonly seen in previous experiments. This will be further discussed in Section 5.1. The probability distribution of the transverse particle velocity, Vy, was plotted in Figure 5. Similarly, Vy shows a symmetric bell-shaped distribution that has a narrow central peak and heavy tails on both sides that are much thicker than those of the Gaussian distribution. Thus, we proposed that the α-stable distribution also provides a better fit for Vy with proper parameters. Our experimental results contradict the widely accepted view that the transverse particle velocity follows the thin-tailed exponential or Gaussian distributions based on previous experiments [10,13,17]. After comparing with the previous experiments, we attribute this inconsistency to the armored gravel bed configuration with cluster microforms, which was not commonly seen in previous experiments. This will be further discussed in Section 5.1.

Particle Step Length and Rest Period
The particle step length, s, was fitted by the Gaussian distribution and the α-stable distribution, which were plotted in Figure 6. Because there are a few negative particle steps in our dataset, the two-sided Gaussian and α-stable distributions were adopted. As shown in Figure 6, the heavy-tailed α-stable distribution with power-law tails fits the experimental data better than does the thin-tailed Gaussian distribution. Nakagawa [46] reported that the distribution of the particle step length was thin-tail for uniform-sized particles. Ganti et al. [26] and Hill et al. [29] suggested that the distribution of the particle step length for mixed-sized sediment was heavy-tailed because of the size heterogeneity in moving grains. In our experiment, the tracer particles were approximately uniform in size; therefore, the heavy tails in the distribution of the particle step length is unlikely to be caused by size heterogeneity. It is noted that the probability distribution of the particle length is similar to

Particle Step Length and Rest Period
The particle step length, s, was fitted by the Gaussian distribution and the α-stable distribution, which were plotted in Figure 6. Because there are a few negative particle steps in our dataset, the two-sided Gaussian and α-stable distributions were adopted. As shown in Figure 6, the heavy-tailed α-stable distribution with power-law tails fits the experimental data better than does the thin-tailed Gaussian distribution. Nakagawa [46] reported that the distribution of the particle step length was thin-tail for uniform-sized particles. Ganti et al. [26] and Hill et al. [29] suggested that the distribution of the particle step length for mixed-sized sediment was heavy-tailed because of the size heterogeneity in moving grains. In our experiment, the tracer particles were approximately uniform in size; therefore, the heavy tails in the distribution of the particle step length is unlikely to be caused by size heterogeneity. It is noted that the probability distribution of the particle length is similar to that of the particle streamwise velocity, which is not surprising because both parameters were related to the preferential particle transport and trapping mechanisms of particle clusters. Therefore, the heavy-tailed distribution of the particle step length is also likely attributed to the bed configuration in our experiment. The distribution of the rest period was fitted by the exponential distribution and the M-L distribution in Figure 7. The M-L distribution with β ∈ [0, 1] is a family of one-sided heavy-tailed probability distributions with power-law tails, which grow thicker as β decreases. Figure 7 shows that the M-L distribution fits the data better than the exponential distribution, which indicates that a portion of the bed load particles experienced unusually long rest periods than expected by the exponential distribution. This is consistent with the results in several laboratory experiments [11] and field surveys [28], which suggested that the heavy-tailed distribution of the particle rest period occurs because particles were buried beneath the surface or stored within the bars and flood plains of a river [47]. In our experiment, particle clusters may obstruct the transport of bed load particles and increase the chances of particle trapping, thereby resulting in usually long particle rests in bed load transport. Figure 7. (a) Discrete probability distribution of the particle rest period (Δt = 0.1 s) with theoretical lines that fit the data; the inset shows the same data on a semilog plot. (b) Exceedance probability of the particle rest period. The distribution of the rest period was fitted by the exponential distribution and the M-L distribution in Figure 7. The M-L distribution with β ∈ [0, 1] is a family of one-sided heavy-tailed probability distributions with power-law tails, which grow thicker as β decreases. Figure 7 shows that the M-L distribution fits the data better than the exponential distribution, which indicates that a portion of the bed load particles experienced unusually long rest periods than expected by the exponential distribution. This is consistent with the results in several laboratory experiments [11] and field surveys [28], which suggested that the heavy-tailed distribution of the particle rest period occurs because particles were buried beneath the surface or stored within the bars and flood plains of a river [47]. In our experiment, particle clusters may obstruct the transport of bed load particles and increase the chances of particle trapping, thereby resulting in usually long particle rests in bed load transport. that the M-L distribution fits the data better than the exponential distribution, which indicates that a portion of the bed load particles experienced unusually long rest periods than expected by the exponential distribution. This is consistent with the results in several laboratory experiments [11] and field surveys [28], which suggested that the heavy-tailed distribution of the particle rest period occurs because particles were buried beneath the surface or stored within the bars and flood plains of a river [47]. In our experiment, particle clusters may obstruct the transport of bed load particles and increase the chances of particle trapping, thereby resulting in usually long particle rests in bed load transport. Figure 7. (a) Discrete probability distribution of the particle rest period (Δt = 0.1 s) with theoretical lines that fit the data; the inset shows the same data on a semilog plot. (b) Exceedance probability of the particle rest period.

Comparison with Previous Experiments
The present experiment has some important differences with previous bed load transport experiments. A distinct feature is the use of armored bed layer with cluster microforms, which was formed under constant flow with widely graded sediment particles. As we reviewed earlier, most of the previous experiments used either uniformly sized [13,17] or narrowly graded sediment [10,12,29]

Comparison with Previous Experiments
The present experiment has some important differences with previous bed load transport experiments. A distinct feature is the use of armored bed layer with cluster microforms, which was formed under constant flow with widely graded sediment particles. As we reviewed earlier, most of the previous experiments used either uniformly sized [13,17] or narrowly graded sediment [10,12,29] as bed material, in which the particle clusters were rarely encountered. Particle clusters may alter bed particle organization and local turbulence structure, and thus has a strong impact on the bed load particle motions. Consequently, the bed load particle motions exhibit some unique features that have not been reported in previous experiments. For example, the probability distributions of both the streamwise and the transverse particle velocities were heavy-tailed, indicating that a considerable amount of particles travelled significantly faster than those predicted by the Gaussian or the exponential distributions proposed in previous studies. The discrepancy is likely a result of the bed configuration in our experiment. Figure 8 shows a picture of the quasi-equilibrium armored bed surface with cluster microforms. The largest particles with a medium size of 17 mm (marked by red circles) were observed to congregate and form cluster structures of different shapes, including the line, heap, and ring clusters. When the tracer particles travelled through the paths between the clusters (marked by blue lines), they were accelerated due to locally enhanced water flow, thereby producing the heavy tails in the streamwise particle velocity distribution. Brayshaw et al. [48] measured the pressure field around the anchor stone on sediment beds and reported that the flow was enhanced on the flanks of the anchor stone to produce two symmetrical high shear stress regions in which fine particles are unlikely to rest. Note that some of the particle paths were highly tortuous due to the obstruction of cluster structures; therefore, the probability distribution of the transverse particle velocity also was also affected by these accelerating particles and became heavy-tailed. clusters (marked by blue lines), they were accelerated due to locally enhanced water flow, thereby producing the heavy tails in the streamwise particle velocity distribution. Brayshaw et al. [48] measured the pressure field around the anchor stone on sediment beds and reported that the flow was enhanced on the flanks of the anchor stone to produce two symmetrical high shear stress regions in which fine particles are unlikely to rest. Note that some of the particle paths were highly tortuous due to the obstruction of cluster structures; therefore, the probability distribution of the transverse particle velocity also was also affected by these accelerating particles and became heavy-tailed. Another notable feature of the bed load particle motions is the high fraction of negative streamwise particle velocity. In previous experiments, negative streamwise particle velocities were treated as "false" velocities due to water surface oscillation [10] or trembling motions of resting particles [17], and were excluded from the data sets, and the rest of the data were fitted by one-sided probability distributions. In our experiment, negative velocities accounted for almost 10% of the total Vx measurements, which is unlikely to be measurement error only. During the experiment, the tracer particles were occasionally observed to quiver inside the ring clusters, or bounce back when colliding with large resting particles in the armored bed layer. Therefore, we proposed that negative streamwise particle velocities should be considered in studying bed load sediment transport.
Based on our observations, we argue that the particle clusters over armored bed layers exert a dual impact on bed load sediment particles: some particles are accelerated in the preferential paths between particle clusters, while others were obstructed by the clusters to produce negative streamwise particle velocities. The particle clusters may also influence the probability distributions of other kinematic parameters. For example, particle acceleration in preferential transport paths between clusters and particle trapping in the ring-shaped clusters may increase the opportunity of Another notable feature of the bed load particle motions is the high fraction of negative streamwise particle velocity. In previous experiments, negative streamwise particle velocities were treated as "false" velocities due to water surface oscillation [10] or trembling motions of resting particles [17], and were excluded from the data sets, and the rest of the data were fitted by one-sided probability distributions. In our experiment, negative velocities accounted for almost 10% of the total Vx measurements, which is unlikely to be measurement error only. During the experiment, the tracer particles were occasionally observed to quiver inside the ring clusters, or bounce back when colliding with large resting particles in the armored bed layer. Therefore, we proposed that negative streamwise particle velocities should be considered in studying bed load sediment transport.
Based on our observations, we argue that the particle clusters over armored bed layers exert a dual impact on bed load sediment particles: some particles are accelerated in the preferential paths between particle clusters, while others were obstructed by the clusters to produce negative streamwise particle velocities. The particle clusters may also influence the probability distributions of other kinematic parameters. For example, particle acceleration in preferential transport paths between clusters and particle trapping in the ring-shaped clusters may increase the opportunity of usually large particle leaps and long particle rests and results in heavy tails in the probability distributions of both the particle step length and the rest period. Nonetheless, the heavy-tailed distributions for the particle step length and the rest period were also reported in several studies where the armored bed layer was absent. More experimental data are warranted to confirm the influence of the armored bed layer with cluster microforms on the heavy tails of these distributions.

Particle Diffusion Regime and Particle Motion Statistics
The scaling growth of the MSD with time indicates that the particle diffusion regime in our experiment is anomalous. In our experiment, the particle diffusion regime was superdiffusive at small time scales, but became subdiffusive at larger time scales. Previous studies suggested that the particle diffusion regime was constrained by the statistical characteristics of particle motions [13,22,38,39]; therefore, the anomalous diffusion regime in our experiment is likely associated with the heavy-tailed distributions of the kinematic parameters. Weeks et al. [36] and Fan et al. [15] investigated the quantitative relationship between the scaling exponent (α) of the MSD and the power-law exponents of the rest period distribution with random walk and Langevin-Fokker Planck Equation, respectively, and their results were rather consistent. To compare with their results, we used the power-law distribution to fit our experimental data, which produced a tail parameter of ν = 1.065 (corresponding to a total of 10 sample groups) for the rest period distribution. According to Weeks et al. [36] and Fan et al. [15], this corresponds to a diffusion exponent of α = 0.13 (strongly subdiffusive) at large time scales, which is very close with the scaling exponent of α = 0.14 computed from the relationship between the MSD and time at the time scales of t > 150 s (Figure 2). It is indicated that, despite a heavy-tailed distribution of the particle step length, the diffusion regime mainly depends on the rest regime rather than the motion regime of the bed load particles.
The particle transport regime may be scale-dependent, i.e., the transport process may exhibit different diffusion regimes on different timescales. The multiplicity of temporal scales in bed load transport has been noted in the literature [15,39,49,50]. Here we discussed the scale dependence of bed load motion parameters using the step length as an example. Figure 9 compares the best-fit Gaussian and α-stable distributions of the step length for observation intervals of ∆t = 0.1, 0.3, 0.5 and 1 s. The shape parameters of the fitted α-stable distributions are 1.290, 1.326, 1.384, and 1.727 for ∆t = 0.1, 0.3, 0.5 and 1 s, which indicates that the tails of the α-stable distributions decay faster and the gradually approaches the tails of the Gaussian distributions as ∆t increases. Similarly, the probability distributions of the streamwise and transverse particle velocities and the particle rest period (not shown) all exhibit thinner tails as observation time scale increases. Voepel et al. [28] suggested that the particle rest period distribution changes from a heavy-tailed distribution to an exponential distribution as the time scale increases. Therefore, we proposed that the heavy-tailed distributions in kinematic parameters in bed load transport are associated with small time scales and may average out over longer observation periods. Nonetheless, we should notice the difference between the time scales associated with the transition of particle diffusion regimes and those with the changes in the statistical characteristics of bed load particle motions. Further investigation is needed in unraveling their connections.
particle diffusion regime was constrained by the statistical characteristics of particle motions [13,22,38,39]; therefore, the anomalous diffusion regime in our experiment is likely associated with the heavy-tailed distributions of the kinematic parameters. Weeks et al. [36] and Fan et al. [15] investigated the quantitative relationship between the scaling exponent (α) of the MSD and the power-law exponents of the rest period distribution with random walk and Langevin-Fokker Planck Equation, respectively, and their results were rather consistent. To compare with their results, we used the power-law distribution to fit our experimental data, which produced a tail parameter of ν = 1.065 (corresponding to a total of 10 sample groups) for the rest period distribution. According to Weeks et al. [36] and Fan et al. [15], this corresponds to a diffusion exponent of α = 0.13 (strongly subdiffusive) at large time scales, which is very close with the scaling exponent of α = 0.14 computed from the relationship between the MSD and time at the time scales of t > 150 s (Figure 2). It is indicated that, despite a heavy-tailed distribution of the particle step length, the diffusion regime mainly depends on the rest regime rather than the motion regime of the bed load particles.
The particle transport regime may be scale-dependent, i.e., the transport process may exhibit different diffusion regimes on different timescales. The multiplicity of temporal scales in bed load transport has been noted in the literature [15,39,49,50]. Here we discussed the scale dependence of bed load motion parameters using the step length as an example. Figure 9 compares the best-fit Gaussian and α-stable distributions of the step length for observation intervals of Δt = 0.1, 0.3, 0.5 and 1 s. The shape parameters of the fitted α-stable distributions are 1.290, 1.326, 1.384, and 1.727 for Δt = 0.1, 0.3, 0.5 and 1 s, which indicates that the tails of the α-stable distributions decay faster and the gradually approaches the tails of the Gaussian distributions as Δt increases. Similarly, the probability distributions of the streamwise and transverse particle velocities and the particle rest period (not shown) all exhibit thinner tails as observation time scale increases. Voepel et al. [28] suggested that the particle rest period distribution changes from a heavy-tailed distribution to an exponential distribution as the time scale increases. Therefore, we proposed that the heavy-tailed distributions in kinematic parameters in bed load transport are associated with small time scales and may average out over longer observation periods. Nonetheless, we should notice the difference between the time scales associated with the transition of particle diffusion regimes and those with the changes in the statistical characteristics of bed load particle motions. Further investigation is needed in unraveling their connections.

Limitations of the Present Study
The present study may suffer from some limitations due to relatively low sampling frequency and a limited observation window. The sampling frequency determines the smallest possible time scale that can be detected in our image analysis and therefore has a significant influence on the

Limitations of the Present Study
The present study may suffer from some limitations due to relatively low sampling frequency and a limited observation window. The sampling frequency determines the smallest possible time scale that can be detected in our image analysis and therefore has a significant influence on the analysis results. Currently, there is no consensus on the requirement of the minimum time scale. Theoretically speaking, the sampling frequency should be high enough to match the smallest step time or rest period of bed load particles; in other words, the time interval between two frames should be equal or smaller than the transition time between two particle states. Our experiments deal with low-intensity bed load transport in which the rest period accounts for 98% of the total observation time. Increasing the existing sampling frequency (10 frames per second) may potentially thicken the distribution tails of particle velocity and step length, but is unlikely to affect the heavy-tailed nature of these distributions.
The mass loss of tracer particles due to a limited observation window size may induce measurement errors that lead to misjudgments about particle diffusive regimes [14,51]. In our experiment, because of the low transport flux and strong resistance of the armored bed layer to entrainment, only about 1% of the tracer particles left the sampling pool during the study period. Hence, the measurement error due to the window size is not significant. Nonetheless, for experiments with large sediment flux, small observation windows may result in false conclusions due to sampling bias towards small step lengths. In the study of particle anomalous diffusion, attention should be paid to the "window" effect on the calculated variance of particle positions. To further explore the details in anomalous bed load transport, higher temporal resolution and larger observation area may be needed in future experiments.

Conclusions
This study aims to explore the stochastic nature of bed load transport over an armored bed with cluster microforms. Distinct from previous experiments in which uniform or narrowly graded particles were used, our experiment employed widely graded sediment particles to form the armored bed layer under constant flow. The statistical properties of a series of kinematic parameters, including the MSD, the particle number activity, the streamwise and transverse particle velocities, the step length, and the rest period, were analyzed, and the diffusion regime of the bed load transport was discussed. The major findings of this study are summarized as follows: 1.
The scaling growth of the MSD with time showed that the particle diffusion regimes in this experiment varied at different time scales: the particle diffusion regime was superdiffusive at small time scales (T < 50 s), but became subdiffusive at larger time scales. 2.
The particle number activity followed a negative binomial distribution, which is consistent with the theoretical solution derived in Ancey et al. [22] by accounting for the collective entrainment processes. The power-law tails in the probability distributions of both the streamwise and transverse particle velocities indicated a significant amount of fast-travelling particles, whereas the negative tail in the streamwise particle velocity distribution indicated a noticeable amount of grains occasionally moved in the upstream direction, which has rarely been considered in previous experiments. This may be explained by the particle clusters over dual impact on bed load transport exerted by the armored bed layers: some particles are accelerated in the preferential paths between particle clusters, while others were obstructed by the clusters to produce negative streamwise particle velocities. The distributions of the transverse particle velocity, the step length and the rest period all displayed heavy tails that were well described by the α-stable and M-L distributions. Nonetheless, the effect of the armored bed layer with cluster microforms on the distributions of these parameters needs to be confirmed with more experimental data.

3.
The particle transport regime may be scale-dependent, i.e., the transport process may exhibit different diffusion regimes on different timescales. Meanwhile, the shapes of the distribution tails of the bed load motion parameters varied as observation time scale increases. We propose that the heavy tails in the distributions of the bed load motion parameters were associated with small time scales and may average out over longer observation periods. Nonetheless, the time scales associated with the transition of particle diffusion regimes and those with the changes in the statistical characteristics of bed load particle motions were different. Further investigation is needed in unraveling their connections.