The transport of sediment mixtures examined with a birth-death model for grain-size fractions

Bedload transport of sediment mixtures in mountain streams is challenging to predict, with implications for understanding how rivers form and respond to environmental change. Experimental work shows that collective particle entrainment is an important contributing mechanism of bedload transport, but questions remain. We use four different time series of experimental sediment flux for granular particles 4–32 mm in diameter to indirectly examine the role of collective mobilization. Flux was measured at a fixed position in space using an imaging light table. The light table provides a flux measurement that 5 is sampled at a resolution of 1 Hz, and for total time durations ranging from 75 to 240 min. Experimental conditions include periods of statistical steady-state, and transient adjustments due to changes of the upstream supply of water and sediment. We find that despite the contrasting experimental conditions, the time series encode a consistent transport behaviour within the Fourier domain: the transport of finer grain size populations has increasing power density for decreasing frequency, whereas the transport of larger grain size populations has a near constant power density across all frequencies. Hence, smaller particle 10 sizes dominate the power spectra. We seek an explanation for this result, and elaborate on a probabilistic birth-death model introduced to the field by Christophe Ancey and colleagues. Analysis using the expanded birth-death model provides two important results. The transport of smaller particles includes collective entrainment terms that represent grain mobilization due to smaller and larger particle sizes colliding with the streambed surface. In contrast, the transport of larger particles includes collective entrainment terms limited to larger particle sizes. The size-dependent collective controls on particle mobilization is 15 an important finding, and we show that it offers a testable explanation for observed flux differences between smaller and larger particle sizes, common to gravel-bed mountain streams. As a result, our work motivates the need to better understand collective entrainment within the context of granular sediment transport along mountain stream beds.

tenet of this work is that the entrainment or mobility of particles resting on a river bed is a stochastic process, with distributions of travel distances and inter-travel period resting times (see Einstein, 1937;Hassan et al., 1991). Consequently, the number of moving particles per unit area of the bed, and the particle motion components of travel time, travel distance and downstream average velocity are characterized by probability distributions with right tails (Hubbell and Sayre, 1964;Hassan et al., 1991;Lajeunesse et al., 2010;Furbish et al., 2012b;Roseberry et al., 2012;Fathel et al., 2015). 15 Experimental observations of particle motions have typically involved tracking grains moving within an image reference area. Here we ask what can be learned from the analogous case of measuring sediment flux at a fixed position in space, at a high temporal resolution and for a period of time that is long relative to particle travel times of fractions or a few seconds or more (e.g. Drake et al., 1988;Roseberry et al., 2012;Fathel et al., 2015). More specifically, what do near-continuous time series of sediment flux reveal about the transport of sediment mixtures? We know of only a limited number of studies that have 20 reported and analyzed near-continuous measurements of sediment flux at a fixed position, and for time scales that range from minutes to hundreds of hours (e.g. Frey et al., 2003;Zimmermann et al., 2008;Singh et al., 2009;Jerolmack and Paola, 2010;Turowski et al., 2011;Ferrer-Boix and Hassan, 2015;Saletti et al., 2015;Dhont and Ancey, 2018;Lee and Jerolmack, 2018;Phillips et al., 2018;Redolfi et al., 2018;Masteller et al., 2019). Dhont and Ancey (2018) analyze a long time series of sediment flux data from a flume outlet. Their experimental design 25 emphasized the examination of how bed morphology affects sediment transport under steady conditions of upstream sediment supply. Their experimental duration was approximately 550 h, during which sediment output from the flume was averaged over periods of 1 min. Sediment flux shows variations that span several orders of magnitude around the sediment supply rate. Dhont and Ancey (2018) report that excursions above the supply rate were typified by pulses which rise to a maximum value over relatively short time intervals. The intermittent sharp increases of sediment flux within the time series results in a 30 visually erratic signal (Singh et al., 2009). However, the flux power spectrum suggests order within the signal, notably that flux variations with frequencies f between 10 −4.8 to 10 −2 Hz fall off as f −3/2 Hz. Dhont and Ancey (2018) indicate that the lower frequency 10 −4.8 Hz flux component corresponds to the approximate arrival of migrating bars, leading to a majority of the measured peak bedload transport rates. The higher frequency flux variations occur within the bar inter-arrival times. The fall off slope of the calculated flux power spectrum reported by Dhont and Ancey (2018) is a striking result of their work. The authors did not provide an explanation for possible controls on the spectra fall off slope. As we discuss below, this behaviour is a primary focus of our work. Saletti et al. (2015) take a different approach and use 13 sediment sampling intervals lasting 10 min to evaluate sediment transport dynamics from an experimental step-pool channel. Sediment flux measured at the flume outlet was the outcome of channel adjustment in the absence of upstream sediment supply, coupled with a stepped hydrograph marked by periods of 5 steady water supply. The authors consider sediment flux as an autocorrelated phenomenon, and show that the Hurst exponent (H) exhibits dependence on the sediment flux magnitude and grain size. The reported experimental conditions of larger flux magnitudes and finer flux compositions have H → 1, whereas smaller flux magnitudes and coarser compositions have H → 0.5. These results imply that the transport of finer sediments upstream of the flux measurement point is more continuous in time and space (also see Kuhnle and Southard, 1988). This contrasts with the transport of coarser particles, which by comparison 10 is more intermittent. Saletti et al. (2015) attribute the different flux responses to the time varying magnitude of upstream bed slope adjustment. Lee and Jerolmack (2018) report results from a series of experiments designed to examine how bedload transport responds to changes in the driving frequency of transport, controlled by different sediment supply rates to a small flume. High-speed imaging of particle motions shows that entrainment of marbles occurs through collisions of moving particles with the bed 15 surface, generally setting two or more marbles into collective motion. The total displacement distance of marbles, and the mean number of marbles set into motion scales in a linear way with the total kinetic energy transferred to surface particles by colliding marbles. However, entrainment of a given number of marbles is shown to be described by a distribution of total kinetic energies, suggesting that entrainment depends on numerous local conditions, which are difficult to anticipate. Related to our work here, Lee and Jerolmack (2018) provide clear evidence that collective entrainment is an important physical process 20 of bedload transport. This finding links directly with theoretical work by Ancey et al. (2008) and Ancey (2010) who develop a stochastic model of bedload transport that explicitly accounts for collective particle entrainment.
Here, we build upon the work of Ancey et al. (2008) and Ancey (2010), as elaborated by  and Heyman et al. (2014), Saletti et al. (2015) and Dhont and Ancey (2018), and examine sediment transport through the information embedded within an experimental particle flux time series of two contrasting conditions: steady-state followed by a transient. In our work, 25 steady-state is defined by two different conditions for periods of time much greater than time averaged particle travel times: (1) the sediment supply rate is comparable to the flux, and (2) changes to the mean bed elevation approaches zero. In contrast, a transient state is characterized by relatively large negative differences between the rate of sediment supply and flux, therefore the channel profile undergoes rapid changes. Our focus is important because the time varying character of sediment flux yields specific information about the upstream bedload transport process. Furthermore, an examination of flux time series between 30 steady-state and transient conditions permits us to explore whether bedload transport processes depend on the adjustment state of upstream areas of the river bed which contribute directly to the measured flux. Furthermore, our experiments include initial and repeat phases, providing for a comparison of flux between specifically two sets of overall similar conditions. We address several questions. First, are all grain size fractions transported at similar rates determined by time averaged counts of moving particles? It is generally accepted that the relative proportion of grain sizes present on the bed surface will influence the overall flux magnitude for each particular size class (e.g. Parker, 1990;Wilcock and McArdell, 1997;Wilcock and Crowe, 2003). However, are effects related to grain size proportion evident in the time series signal of grain size specific flux?
Second, do statistical steady-state and transient sediment signals carry the same information about the transport process? That is, other than possible differences in flux magnitude, can we specifically infer from flux data if the bedload transport process is similar, or differs between steady-state and transient? Third, the transport process information encoded within sediment flux 5 time series can be understood through power spectral analysis. Jerolmack and Paola (2010) and Dhont and Ancey (2018) show that the power spectra of total flux can exhibit, for example, power law structure over a finite frequency band or time scale range. What transport mechanisms are responsible for the spectral structure?
We use time series analysis to explore these questions, and we expand the Markov birth-death framework of Ancey et al. (2008) and Ancey (2010) to particle size mixtures. Our work provides further evidence that sediment mixtures exhibit size-10 dependent transport at low to moderate rates of particle flux. Power spectra for steady-state and transient conditions indicate that fine grain flux is dominated by frequencies less than approximately f −2 Hz. In contrast, coarse particle flux shows no dominant frequency or frequency band, reflecting a white noise bedload transport process. We show that the fine grain flux spectra likely reflects the effects of collective entrainment in the birth-death model, which accounts for larger grains collectively entraining finer ones. This collective effect diminishes with grain size as the pool of particles capable of setting larger grains into motion 15 is presumed to decrease. Our combined results further motivates the critical need for continued development of a probabilistic theory of sediment transport for gravel-bed rivers. Fruitful next steps include an emphasis on further examination of the physical controls of collective particle entrainment. The recent work by Lee and Jerolmack (2018) Figure 2. Definition diagram of a light table system within which the sediment particle flux is measured downstream of a flume (see Zimmermann et al., 2008;Chartrand et al., 2018). The vertical dashed gray lines over the light table mark 5 sub-regions of the light table within which particle flux is estimated using successive images of the light table (15)(16)(17)(18)(19)(20).

Problem Set-up
Our work explores the time dependent character of grain-size specific sediment flux measured in the laboratory during development and adjustment of pool-riffle like bed topography (see Chartrand et al., 2018). Since sediment flux is measured in a variety of ways in the laboratory (e.g. Frey et al., 2003;Zimmermann et al., 2008;Singh et al., 2009) and the field (e.g. Hubbell, 1964;Helley and Smith, 1971;Diplas et al., 2008;Rickenmann, 2018), it is important to establish context for the information 5 carried by the transport time series that we report here (cf. Heyman et al., 2013).

Particle Activity
Consider an idealized streambed area A with streamwise length L and width B. Particles enter A through the left boundary and exit through the right boundary. Neither entrainment nor deposition occurs within A. At any time t there are N (t) moving particles within A and the associated particle number activity is γ(t) = N (t)/A. 10 Let n(t) denote the number of particles of a given size that exits the right boundary following an arbitrary starting time t = 0 ( Fig. 1). The number discharge specific to a measurement interval ∆t is then where for convenience we associate this finite difference with time t. Because n(t) is a stepped function of time, the limit of the right side of Eq. (1) as ∆t → 0 is strictly zero or undefined. Thus Q ∆t represents an average discharge over the interval ∆t. 15 Consider the light table configuration described by Zimmermann et al. (2008) (Fig. 2 If between time t and time t + ∆t a number of particles, ∆n, crossed the downstream boundary of the light table, then assuming a constant particle velocity u * p over the table, this number of particles resided on the light table within some area A * at time t. At this instant the particle activity associated with A * is γ * (t) = ∆n/A * . That is, upon setting ∆n = N (t), this is 5 the activity γ * (t) = N (t)/A * on the light table, although A * may differ from the actual area of the apparatus.
The vertically integrated flux is q ∆t (t) = Q ∆t (t)/B * . This implies that q ∆t (t) = γ * (t)L * /∆t. Setting u * p = L * /∆t then gives where L * now must be interpreted as the distance traveled by particles over the light table at velocity u * p during the interval ∆t.

10
The number of particles ∆n that cross the right boundary of the light table during ∆t and which resided on the light table at time t were delivered from upstream across the left boundary of the light table (Fig. 2). By continuity, this influx across the left boundary is equal to that expressed in Eq.
(2), but it occurred during an interval ∆t prior to time t. By an argument similar to that above, we can envision a streambed area A contributing to this flux equal to 15 where γ is the number activity associated with A and u p is the average particle velocity. If variations in the velocities of a given size are small, then the distance L is interpreted as the distance traveled by particles over the streambed at velocity u p during ∆t. Moreover, if we momentarily assume that the difference between entrainment and deposition within A during ∆t is negligible, then 20 That is, the particle number activity within the streambed area A upstream from the light table is directly proportional to that measured on the light table, albeit separated in time by at least ∆t.
In fact, entrainment and deposition within an area A are not likely to be balanced during ∆t. What is important is that the activity γ * is a measure of the upstream activity in Eq.
(3) of particles delivered to the light table. Our association of upstream activity and flux measured by the light table is an important assumption of our expanded birth-death formulation, which we 25 now present.

Birth-Death Formulation
Following Ancey (2010), for a specified streambed area A (Fig. 2), the number of active particles of the jth size class may be expressed as Here, v j is the immigration rate, λ j is the entrainment rate normally attributed to fluid forces, µ j is a rate constant of collective entrainment, σ j is a rate constant of deposition, and w j is a rate constant of emigration. To be clear, the fluid related entrainment rate is viewed as dislodgement from the bed due to fluid forces, as in the classic assumption in sediment transport (e.g. Shields, 1936;Meyer-Peter and Muller, 1948). Collective entrainment is viewed as dislodgement "because of the moving particles, which can interact directly (collision) or indirectly (wake effect, advection of turbulent structure) with the bed particles" 5 (Ancey, 2010). To this we add the possibility that entrainment of a particle from the bed leads to the onset of motion of nearby particles with localized rearrangement and destabilization of the bed (Böhm et al., 2004;Lee and Jerolmack, 2018;Chartrand, 2017).
Letting a circumflex denote a Fourier transform, the birth-death formulation, Eq. (5), is represented in the frequency domain as: where ω = 2π/T is the angular frequency with period T . For a specified area A, the immigration rate v j and the emigration rate w j n j may differ at any instant. However, because the emigration rate for any area A is just the immigration rate for the downstream area, the spectra of these must be the same for statistically uniform transport. This means thatv j = w jnj , and we rewrite Eq. (6) (see Appendix D): Rearranging Eq. (7) provides the grain size specific Fourier transform of the number of active particles for a specified streambed area A (Fig. 2): The Fourier transformed particle activityn j is determined by the ratio of two different physical affects. First,n j depends 20 directly on the fluid entrainment rates across the frequencies iω. For example, at low to moderate bedload transport stages (Hassan et al., 2005), particle entrainment frequencies due to the time-averaged downstream fluid velocity versus that due to turbulent sweeps of fluid at the bed differs. The former produces fluid entrainment of particles at frequencies that are high relative to turbulent sweeps, which occur at some characteristic time scale (Ancey et al., 2008;Singh et al., 2009).
Second,n j depends inversely on the balance between particle deposition and collective entrainment rate constants, which can 25 attenuate or amplifyλ j . In general,λ j is amplified at times when deposition and collective entrainment are comparable, and is attenuated when deposition is large relative to collective entrainment, assuming µ j σ j (e.g. Ancey, 2010;Heyman et al., 2014). However, the amplification or attenuation effect toλ j is influenced further by the range of frequencies iω of bedload transport in gravel-bed rivers. As a result, regardless of whetherλ j is amplified or attenuated by (σ j − µ j ) −1 , the sum acts as a low pass filter onλ j . Physically, this result indicates that bedload transport is characterized by relatively high frequency 30 fluctuations (e.g. Ancey et al., 2006;Jerolmack and Paola, 2010;Fathel et al., 2015). We discuss these issues in more detail below following further development of the birth-death model.
Log frequency ( Figure 3. Example of how an input white noise for n3 contributes to a red noise signal for n1 based on our imagined 3 particle size system, and over the frequency band 10 −3 to 10 −0.5 Hz. The power spectral density of (a) is based on randomly sampling a Poisson distribution Expanding on Eq. (5), consider three particle size classes (j = 1, 2, 3), small to large. For simplicity, assume that collective entrainment of the jth size involves the jth and larger sizes: dn 2 dt = λ 2 + µ 22 n 2 + µ 32 n 3 − σ 2 n 2 , and 5 dn 1 dt = λ 1 + µ 11 n 1 + µ 21 n 2 + µ 31 n 3 − σ 1 n 1 .
The first subscript on the collective entrainment rate constant µ denotes the influencing particle size and the second subscript denotes the influenced size. Taking the Fourier transforms, Rearranging then leads to: Equations 15-17 indicate that the spectrum of the jth size involves the sum of the spectra of coarser sizes, each modified by a low-pass filter of the same form, i.e. 1/[σ j − µ ii + iω]. We simplify the notation and use F j to represent the low pass filters: 20 In order to solve Eqs. (18)(19)(20) and learn about the three particle size birth-death model, we need to complete two last steps.
First, we find the real part of the low pass filters Re(F j ) (see Appendix D): Log frequency (

Re(F
Next, we need to specify or calculate the Fourier transform of the time varying fluid entrainment signals for each particle size 5λ j , Eq. (5) and Eqs. (18)(19)(20). Based on experimental evidence, we assume the largest grain size classλ 3 is a white noise (Saletti et al., 2015). On the other hand, transport of the two smaller grain size classesλ 2 andλ 1 is assumed to follow an autoregressive process, and specifically a second-order AR (2) process. An implication of our assumption is that the fluid entrainment rate of particles at time t is correlated to rates during the two previous times, and furthermore may carry longer term memory of prior particle entrainment events (Ancey et al., 2008;Saletti et al., 2015). The temporal correlation is likely stronger, however, for 10 all but the weakest or most sporadic transport conditions (e.g. Drake et al., 1988;Ancey et al., 2008;Roseberry et al., 2012;Fathel et al., 2015;Saletti et al., 2015;Lee and Jerolmack, 2018). As a result, use of an AR(2) process to modelλ 2 andλ 1 is at best illustrative of actual fluid entrainment conditions. The power spectra forλ 2 andλ 1 of an AR(2) process is defined as (Box et al., 2008): The term s 2 j is the variance of an input white noise process, φ aj and φ bj are the autoregressive model coefficients, which here have positive values due to the low-pass behaviour of Eq. (8). Equations (18-25b) along with values for φ aj , φ bj and µ jj permit calculation of the flux power spectral densitiesn 1 -n 3 for our exploratory three particle size system. To that end, we provide two examples to illustrate the behaviour of the birth-death model presented above.

5
Based on Saletti et al. (2015), suppose thatλ 3 is a white noise measured at the upstream boundary of A (Figs. 2 and 3a).
This signal is passed through the low-pass filter F 3 (ω) (Eq. 18), resulting in pink-like noise structure forn 3 (Fig. 3b). Next, n 3 is modulated by µ 32 (Eq. 19) and added toλ 2 . The resulting quantity is then passed through the low-pass filter F 2 (ω) (Eq. 19), resulting in pink noise structure forn 2 (Fig. 3c). The process continues and expands for the finest grain size (Eq. 20), and results in red noise structure forn 1 (Fig. 3d). So, even ifλ 3 is initially a white noise (Fig. 3a), the outcome for smaller grain 10 size classes is a spectrum dominated by low frequencies (Figs. 3b-3d). In fact, because the same filtering process is occurring upstream, the signalsn 1 orn 2 are unlikely to be white noise signals (Fig. 3). In this example we use unique values for each of the parameters of Eqs. (18)-(23), which raises the question how does the birth-death model for a mixture of three particle sizes behave under a broader set of parameter conditions? Power spectra dominated by low frequencies for the two smaller particle sizes is achieved in at least two different ways 15 using the birth-death model with a broader set of parameter values (Fig. 4). First, larger AR(2) weight parameters for the two smaller particle sizes amplifies and shifts power to increasingly lower frequencies as particle size decreases (compare Figs. 4a and 4b). A slightly more pronounced effect is achieved by further increasing collective entrainment for the smaller particle sizes relative to the largest size (compare Figs. 4b and 4c). Finally, amplification of the lower frequencies is reduced, but not eliminated, when the AR(2) weight parameter values are decreased, revealing that collective entrainment can amplify lower 20 frequency spectral power (compare Figs. 4a, 4c and 4d) (cf. Ancey and Heyman, 2014). These combined results suggest that lower frequency spectral power is amplified by enhanced effects due to local fluid driven entrainment, as well as collective entrainment (Ancey et al., 2008;Lee and Jerolmack, 2018). However, the amplifying effect of fluid entrainment may be more important. Regardless, more work is needed to better understand the numerical characteristics of terms in Eqs. (18-23) .

Methods
Experimental data reported here corresponds to pool-riffle experiment 1 (PRE1) conducted at the BioGeoMorphic eXperimental Laboratory (BGMX Lab) at the University of British Columbia, Vancouver, Canada. The experimental flume measures 16 m in total length and 1 m in width, and we use 15 m to conduct the experiments. Water recirculates through the flume via a pump, but sediment does not (Fig. 2). We introduce water to the upstream channel boundary through a series of stacked 5 cm 5 diameter plastic pipes, collectively called the flow normalizer. The normalizer is 1 m in length, or roughly 2w in length, where w is the average channel width, and we use it to establish an initially uniform flow. We introduce sediment to the flume via a speed-controlled conveyor, which dumps particles into a mixing chamber we call the randomizer. The randomizer consists of  The gray shaded region highlights the grain size classes reported here.
a vertical shaft with alternating cross-bars spanning the width of the shaft. As particles fall through the mixing chamber, their pathways are interrupted by the cross-bars, which flings the particles along random trajectories, providing a spatially-random distribution of sediment fall points on the inlet flume bed. The randomizer action provides a spatially and temporally uniform inlet boundary condition, which did not require manual adjustment during the experiments due to development of a sediment pile at the flume inlet.

5
The flume outlet elevation is fixed, and the downstream-most 1.0 m of the flume consists of straight channel walls. We chose this outlet configuration to provide controlled conditions for water and sediment leaving the flume, which pass through the particle imaging light box for flux measurement (Fig. 5a, and discussed in more detail below). Figures 5a and 5b show that the experimental channel consists of downstream varying channel width, which for simplicity, reflects the interpolated width condition between inflection points along a stream segment. We achieve the experimental width conditions by constructing a 10 channel inside the flume with rough-faced veneer-grade D plywood, which has a surface roughness that varies from 1-4 mm, or roughly 0.15 to 0.60 times the geometric mean grain size of the experimental grain size mixture. The grain size distribution of the experiments ranges from 0.5-32 mm (Fig. 6), with a geometric mean size of 7.3 mm (D g ), a D 90 of 21.3 mm, and a geometric standard deviation of 2.5 (σ g ). We used experimental water supply rates Q w of 42, 60 and 80 liters per second (l s −1 ) with associated sediment supply rates of 0.5, 0.80 and 1.0 kilograms per minute (kg min −1 ), respectively. Additional 15 details about the experimental design are provided by Chartrand et al. (2018), and we note that the experimental grain size distribution and water supply rates are similar to those reported by (Kuhnle and Southard, 1988).

Experimental procedure
We start PRE1 from a smoothed-bed, uniform slope condition. Prior to smoothing, the full thickness of sediment in the flume was mixed to establish a random size distribution, and to remove textural heterogeneity related to previous experiments. Figure   7a shows the experimental water supply rate in l s −1 versus time in minutes (a), the sediment supply rate in kg min −1 and the 5-min moving average sediment flux in kg min −1 (b), and the average flume bed elevation (c). PRE1 had a total duration 5 of 79.8 h, which consists of an initial and repeat phase (Fig. 7b). The initial phase extends from t e = 0 min, where t e is elapse time, to 2630 min, and the repeat phase extends from t e = 2630 min, to 4790 min (Fig. 7b). Flow and sediment supply continue at constant values until the total sediment flux approximates the sediment supply rate, and in all cases the fractional flux was comparable to the fractional supply (Table 2; Fig. 7b). The specific water and sediment supply rates for PRE1 are provided in Table 1. We use a ramping up and down period of 4-5 min each time the water supply is raised and lowered to and from the experimental flows of 42, 60 and 80 l s −1 . The repeat phase (Fig. 7b) began from the prevailing channel topographic 5 and bed surface sediment sorting conditions established by the end of the initial phase (Fig. 7b).

Bedload flux
We use a light table to measure bedload flux and ensure mass conservation (Frey et al., 2003;Zimmermann et al., 2008). The light table system uses an overhead camera to measure particle positions in a water column 2 to 3 cm thick. The particles and water pass over a positively-sloping semi-transparent lexan base, which is backlit by a constant-voltage LED panel light 10 measuring 610 mm 2 . Images of the silhouetted particles are captured at 15-20 Hz with an Allied Vision Technology GX2300 CCD camera. The camera uses a Kowa Optimed 16 mm 4/3" megapixel LM16XC lens, which was selected specifically for the GX2300 sensor resolution and imaging distance of the setup. Images are processed with LabView T M code to compute the time-averaged flux for all grain size classes >2 mm at a temporal resolution of 1 Hz (Zimmermann et al., 2008). Here we report fractional flux results for grain size classes >4 mm in units of particles per second (pa s −1 In the following subsection we describe analysis of the light table sediment flux time series for specified time durations of PRE1. The light table represents a fixed location in space, situated at the downstream end of the flume (Fig. 5). As a result, the flux time series reflects sediment transport conditions, and the associated particle number activity which occurs over a finite length at the downstream end of the flume (Eq. (4)). Approximating this length scale is difficult because it varies with the macroscopic fluid conditions (e.g. Lajeunesse et al., 2010). However, we speculate a correlation between light table flux 5 information and particle activity over 1-2 average channel widths of flume length, situated upstream of the flume exit.

Time Series Analysis
We examine time series min (herein referred to as t4070), respectively (Fig. 7b). The corresponding transient flux time series for the initial and repeat experimental phases are from elapsed times 2150-2225 (herein referred to as t2150) and 4310-4385 min (herein referred to as t4310), respectively (Fig. 7b). 20 We use several different methods to learn about and characterize the flux time series. We plot the fractional flux as histograms, calculate Poisson and negative binomial probability mass function (pmf) fits, and calculate descriptive statistics. The Poisson and negative binomial distributions represent steady-state conditions for the particle number activity γ of a streambed area A as t → ∞, depending on whether µ = 0, or µ > 0, respectively (Ancey et al., 2008). Recall our assumption that the particle number activity γ measured within the streambed area A upstream from the light table is directly proportional to γ * measured 25 on the light table (Eq. 4). Consequently, we expect the probability mass distribution of γ * measured on the light table to approximate a Poisson or negative binomial, again depending on whether µ = 0, or µ > 0, respectively. As a result, we use the PRE1 flux to directly test Ancey and colleagues steady-state solution for their birth-death Markov process model, which we expand here for sediment mixtures. We note, however, that the PRE1 flux data used for this test represents both steady-state and transient experimental conditions (Fig. 7). For the transient conditions, we assume that the relatively large sample size 30 N for each grain size class (see Table 3) justifies approximating the binomial distribution of the general birth-death model solution (Ancey et al., 2008) as a Poisson distribution, and as a negative binomial distribution depending on contributions from collective entrainment. We calculate Poison distribution fits using the optimize and curve fit tool boxes of the SciPy Python Library (Jones et al., 2001). Bin middle values and associated histogram results for the experimental flux time series are inputs for the distribution parameters k and λ, respectively. We calculate negative binomial distribution fits using the methods of moments introduced by Pafnuty Chebyshev in 1887 (Goswami et al., 2019). Notably, we calculate p as E[X]/V [X] in order to be consistent with Ancey et al. (2008), who represent p as the probability of observing any particular value γ within a control volume above some bed area A, and at any instant in time. See Chartrand et al. (2020) for the source codes within which the calculations are made.

5
We apply discrete Fourier transform to the fractional flux time series to examine the distribution of spectral density across frequencies from 10 −3 to 10 −0.5 Hz. Each fractional flux time series was preprocessed according to several steps before calculating the power spectral density. Preprocessing began by removing anomalous signals related to memory caching problems.
Anomalous signals are identified with an internal error code related to the circumstance of insufficient image information during calculation of the 1 Hz average flux values for each grain size class. Anomalous signals have 1 Hz average flux values that are 2-3 orders of magnitude higher than nearby rates, and generally persist for 1 to 2 s. A total of 2.47% of the raw flux data for t1910, were anomalous, 1.60% for t4070, 0.89% for t2150 and 0.42% for t4310. After removal, the anomalous values were replaced by the mode of the preceding 20 s time interval, or shorter if the anomalous values occurred at the beginning of a data collection period. We selected 20 s for mode calculation because the longest string of anomalous values was 19 s. All signal processing discussed next was completed using the signal toolbox of the SciPy Python library (Jones et al., 2001). for each of the six different grain size classes using the windowed and detrended data.

30
The light table has uncertainty related to the ability to resolve grain diameter within the backlit photographic frame. This challenge arises from a variety of sources, including the camera/lens resolution, spatial consistency of the light table output lumens and particle edge resolution (see Zimmermann et al., 2008, for other considerations). As a result and based on preexperiment calibration tests, we do not report flux results for grain sizes < 4 mm. The light table images were processed by setting threshold gray scale values to isolate sediment particles against the backlit lexan sheet. For any one experimental interval, one or more image processing runs were needed to achieve as close a match as practical to the total mass measured by the sediment catch-basket. An exact match between light table total mass and catch basket total mass is impossible. For example, the light table mass for each grain size fraction is calculated by multiplying the 1 Hz grains-size specific particle activity by the average mass for each grain size class (Chartrand et al., 2018).

5
Despite these challenges, the light table results show good agreement with the basket captured total and fractional masses ( Table 2). The top row of Table 2 shows that the ratio between light table total mass and catch basket total mass ranged from 0.91-1.07. Rows 2-7 show the ratios of the fractional masses between light table and a sieved sample from the catch basket total mass. Absolute differences between the fractional masses ranges from 0-0.04. These basic descriptions set the limit of our ability to characterize measurement uncertainty of the light table. This is primarily because we use different gray scale 10 threshold values between measurement time intervals when processing the light table images. Threshold values differences were necessary due to changes in the ambient lighting condition of the laboratory.

Approximate Steady-State Flux
Statistical steady-state flux is characterized by a time-varying signal across all grain size classes, and ranges from approximately 0-40 pa s −1 (Table 3; Fig. 8). Flux magnitude and signal continuity are correlated with grain size (Fig. 8). Finer grain size fractions from 4-11.3 mm have particle counts that range from approximately 0-40 pa s −1 (Table 3), and are discontinuous in time. The time discontinuity of the 1 Hz flux increases with particle size.
Flux modes for the finer grain size fractions is zero, and the median is nonzero in only two instances across both time series.
The mean flux, on the other hand, is nonzero for all finer size classes, and the variance extends over a wide range of values (Table 3). However, the mean flux magnitude for the finer size classes is sensitive to the aggregation time interval (Fig. A1). For example, the t4070 mean flux of the 4 mm particle size varies by up to a factor 4 for aggregation time intervals that range from t/20-t, where t is the total time duration of the measurement interval. For t4070 t = 240 min.
Flux magnitude for coarser grain size fractions from 16 to 32 mm, on the other hand, have a range from approximately 0-10 pa s −1 , and have a greater proportion of periods of zero flux (Table 3; Figs. 8 and 9). This last point is most clear in the record for the coarsest grain size fraction (22.6-32 mm), which has periods of zero flux that range from 10 to 100 s or more in 5 duration (Fig. 8). The mode and median of the coarser grain size fractions is 0 (Table 3). However, the mean is nonzero in all cases. The positive skew of the steady-state fractional fluxes decreases with grain size (Fig. 9).
Maximum flux for the finer grain size fractions are associated with short durations of elevated transport that are embedded within an increasing background flux (Fig. 8). During the sharp rises in transport, particle counts increase by a factor 2-3 above the preceding particle counts, and persist at these elevated conditions for seconds to minutes (Fig. 8). This is the case for 10 coarser grain size fractions as well, but the magnitude of the sharp flux increases is smaller than that for finer size fractions by a factor of 2 to 3.
Probability mass function fits to the flux histograms is grain size dependent (Fig. 9). The 4-8 mm grain size classes show better overall agreement with the negative binomial pmf, whereas the 11.3-22.6 mm size classes are in general fit equally well by the Poisson pmf and negative binomial pmf.

Transient Flux
Transient flux is elevated compared to statistical steady conditions and is characterized by a time-varying signal across all grain size classes that ranges from approximately 0-120 pa s −1 (Table 3; Fig. 10). Flux magnitude and signal continuity is correlated with grain size (Fig. 10). Finer grain size fractions from 4-11.3 mm have particle counts that range from approximately 0-120 pa s −1 (Table 3), and flux is discontinuous in time, but less so than the steady-state case. Signal discontinuity increases with 20 grain size.
The flux modes for the finer grain size fractions ranges from 0-5, and the median ranges from 2-17 (Table 3). The mean flux is nonzero for all finer size classes, and the variance extends over a wide range of values (Table 3). However, the mean flux magnitude for the finer size classes is sensitive to the aggregation time interval (Fig. A1). For example, the t4310 mean flux of the 4 mm particle size varies by up to a factor 1.7 for aggregation time intervals that range from t/20-t. For t4310 t = 75 min.

25
Flux magnitude for coarser grain size fractions from 16 to 32 mm, on the other hand, have a range from approximately 0-10 pa s −1 , and have a greater proportion of periods of zero flux (Table 3; Figs. 10 and 11). This last point is most clear in the record for the coarsest grain size fraction (22.6-32 mm), which has periods of zero flux that range from 1 to 10 s or more in duration. The mode and median of the coarser grain size fractions is 0 ( Table 3). The mean flux, however, is nonzero for the coarser size classes. The variance is comparable to the mean (Table 3). Similar to the steady-state results, the positive skew of 30 the transient fractional fluxes decreases with grain size (Fig. 11).
Similar to the steady-sate case, maximum flux for the finer grain size fractions are associated with short durations of elevated transport that are embedded within an increasing background rate of transport (Fig. 10). During the sharp rises in transport, particle counts increase by a factor 2-3 above the immediately preceding particle counts, and persist at these elevated conditions for seconds to minutes (Figs. 10). This is the case for coarser grain size fractions as well, but the magnitude of the sharp flux increases is lower than that for finer size fractions by factors of 2 to 20.
Probability mass function fits to the flux histograms is grain size dependent, and more strongly so compared to the steadystate case (Fig. 11). The 4-8 mm grain size classes show remarkable agreement with the negative binomial pmf, and the fit to examine how the specific grain size spectra compare to the power of the smallest size class. The smallest grain size class 4 mm has the largest magnitude NPSD, and the largest grain size class 22.6 mm has the smallest magnitude NPSD. In between these two bounding grain size classes, the NPSD decreases as grain size increases. Relative to the maximum power of 4 mm, the NPSD span roughly 5 orders of magnitude.
The distribution of spectral density across the range of frequencies also differs according to grain size class. The NPSD for 5 grain size classes 4-11.3 mm have corner frequencies f c of approximately 10 −2.5 Hz. The NPSD decreases proportionally as 1/f − 6 3 to 1/f − 4 3 for frequencies in the approximate range 10 −3 to 10 −1.8 . In contrast, the NPSD for grain size class 16 mm decreases proportionally as 1/f −1 to 1/f − 1 3 , and the NPSD for the 22.6 mm size class is distributed with roughly equal spectral density across all frequencies.
The approximate steady-state distribution of spectral density for the grain sizes evaluated here changes if the input data are 10 represented differently. For example, Fig. B1 shows the NPSD calculated using 1 Hz particle mass estimates. The mass-based normalized spectral density magnitude is size dependent, similar to the particle count based NPSD results (Fig. 12). However, the larger grain size classes account for a greater proportion of the total and maximum spectral density, compared to the smaller size classes, because particle mass is proportional to the diameter cubed. 15 Similar to the steady-state results, the transient NPSD estimates differ according to grain size class (Fig. 13). The smallest grain size class 4 mm has the largest magnitude NPSD, and the largest grain size class 22.6 mm has the smallest magnitude NPSD.

Transient power spectral density
In between these two bounding grain size classes, the NPSD decreases as grain size increases. Relative to the maximum power of 4 mm, the NPSD span roughly 5 orders of magnitude.
The distribution of spectral density across the range of frequencies also differs by grain size class. The NPSD for grain is distributed with roughly equal power across all frequencies. Last, the distribution of transient spectral density changes if particle counts are represented by estimated equivalent masses (Fig. B2). Similar to the steady-state case, the mass-based size dependent distribution of NPSD for transient conditions shows that the larger grain size classes account for a greater proportion of the total and maximum spectral density, compared to the smaller size classes.

5
Autocorrelation of sediment flux varies with grain size class for steady-state and transient conditions (Figs. 14 and 15). Smaller grain sizes show stronger autocorrelation with roughly exponential behavior through lags of ≈250-500 s. Autocorrelation strength and time decreases with grain size, and in general, flux of the largest two grain size classes is uncorrelated in time, and exhibits a white noise type structure. Comparison of the two sets of results reveals that the transient autocorrelation is stronger than the steady-state case for the 4-11.3 mm grain size classes (Figs. 14 and 15).

Transport behaviour of gravel-rich sediment mixtures
Flux measurements for PRE1 made with a light table device indicate bedload transport is characterized by a fluctuating time signal, regardless of grain size (e.g. Iseya and Ikeda, 1987;Kuhnle and Southard, 1988;Frey et al., 2003;Saletti et al., 2015).
However, the grain size specific fluctuations represent different transport behaviours. Grain sizes between 4 and 16 mm exhibit 15 stronger time correlated transport behaviour, whereas grain sizes between 16 and 32 mm exhibit increasingly weaker behaviour . The contrasting transport behaviour for small versus large grains describes both steady-state and transient bed topography conditions (Fig. 7c) (see Chartrand et al., 2018Chartrand et al., , 2019. Physically, a congruence of particle flux behaviours between differing rates of bed topography adjustment and bedload transport stages suggests a common set of physical mechanisms are responsible for bedload transport. Saletti et al. (2015) report similar results, and find that smaller particle sizes tended to exhibit time correlated transport behaviour, whereas increasingly larger grain sizes tended to time uncorrelated behaviour (also see Kuhnle and Southard, 1988).
Consistency between our results and the work of Saletti et al. (2015), and to some extent Kuhnle and Southard (1988), is notable because the experimental conditions differed. The bedload transport data summarized here reflects experiments 5 conducted to examine the development of pool-riffle channel segments within a variable width flume (cf. Chartrand et al., 2018Chartrand et al., , 2019. In contrast, Saletti et al. (2015) report on bedload transport measured during the experimental development of step-pool channel segments, and Kuhnle and Southard (1988) report experiments conducted within uniform width flumes designed to examine bedload transport behaviour. The PRE1 experiments were conducted with steady supplies of sediment at three different water supply rates. Kuhnle and Southard (1988) also used steady upstream supplies of sediment and water, 10 however their sediment supply rates were up to an order of magnitude greater than rates for PRE1. Saletti et al. (2015), on the other hand, report results of experiments for which the water and sediment supply were varied in a more complicated manner.
Despite these differences of experimental design, a consistent size dependent transport behaviour for sediment mixtures is evident. This finding raise an important question: what physical mechanisms might contribute to similar grain size dependent transport behaviour across differing bed architectures, channel slopes and driving conditions? The exploratory three particle 15 size birth-death model generally supported by our experimental results suggests that fluid and collective entrainment, and the interplay between these mechanisms, offers one possible explanation (e.g. Eqs. [15][16][17] (Figs. 3-4, 9 and 11). We elaborate this proposal below by discussing our visual observations of how bedload transport occurred within the flume during PRE1.
We observed several different transport processes during PRE1 that are consistent with previous work, and which together provide context for the birth-death model presented in Section 2. During the steady-state periods t1910 and t4070, particles 20 commonly moved along roughly straight, downstream-oriented transport pathways, including up and to the flume outlet. These pathways were located within the central part of the flume, and the pathways commonly extended for distances of several meters in length (e.g. Powell, 1998;Dhont and Ancey, 2018). The specific location of the pathways might have been controlled by local topographic conditions (e.g. Powell, 1998), but we are not able to demonstrate this effect. Particle sizes up to approximately 16 mm moved within these zones, which often exhibited a consistent transport of grains. As the transport intensity increased 25 during the transient conditions t2150 and t4310, the number of transport pathways increased, and some merged (cf. Lee and Jerolmack, 2018). However, in both cases, it was uncommon to observe more than the occasional grain from the coarsest size classes moving within the pathways. By comparison, the coarsest grains tended to move along pathways that were related to the particles original position relative to the flume walls. We interpret the near time-continuous transport of grains along preferential pathways as a key contribution to the average particle flux for grain sizes up to approximately 16 mm, which is 30 sensitive to the averaging time scale (Figs. A1 and A2) (Kuhnle and Southard, 1988;Bunte et al., 2004;Singh et al., 2009;Saletti et al., 2015;Dhont and Ancey, 2018;Ancey and Pascal, 2020).
We also observed gravel-dune type features many D 50 grain diameters tall (Carling, 1999;Chartrand et al., 2017, 2018, t2150 of the latter), and to a lesser degree bedload sheets a few D 50 grain diameters tall (Whiting et al., 1988;Kuhnle and Southard, 1988;Dhont and Ancey, 2018, e.g.). We did not observe the systematic longitudinal grain sorting reported 35 27 https://doi.org /10.5194/esurf-2021-16 Preprint. Discussion started: 9 March 2021 c Author(s) 2021. CC BY 4.0 License. by Iseya and Ikeda (1987), but the experimental beds were commonly congested in appearance. Both bedload features were observed during steady-state and transient periods, however it was easier to observe the gravel dunes move to the flume outlet intermittently, and result in elevated particle flux that lasted for durations of up to minutes. Unlike the gravel bars reported by Dhont and Ancey (2018), neither gravel-rich feature observed here were consistent attributes of the bed architecture along the flume, which suggests other transport mechanisms such as collective entrainment contributed to the fractional flux measured 5 by the light table (Figs. 8 and 10). We note, though, that the duration of our experiments was short compared to those reported by Dhont and Ancey (2018), and as a result our interpretation of flux-related affects for gravel-dune type features and bedload sheets may be biased.
Based on limitations of our experimental measurements, we cannot assign a specific entrainment mechanism to the downstream motion of sediment particles (Singh et al., 2009). However, the combination of fluid and collective entrainment, and the 10 grain size dependent interplay between these mechanisms provides a plausible explanation. General consistency between the theoretical predictions of the Ancey and colleagues birth-death model (Ancey et al., 2008;Ancey, 2010, andalso see Heyman et al. (2013) and Heyman (2014)), illustrated behaviour of the exploratory three particle size birth-death model developed here , and our experimental results support this perspective. In particular, negative binomial distribution fits to flux histograms for smaller particles size classes of PRE1 is consistent with birth-death formulation involving collective entrainment 15 (Figs. 9 and 11). Going forward, substantial work remains to more explicitly connect mechanisms of entrainment, etc. within gravel-bed river segments to bedload transport processes and associated measurements. To this end, in the next section we discuss testable links between transport fluctuations and the physics of bedload transport, and in our concluding remarks we suggest one path forward to examine identified links.
5.2 What controls the power spectral density of the particle flux? 20 We hypothesize that the power spectra of bedload transport within gravel-bed rivers at low to moderate transport stages (Hassan et al., 2005) can be conceptually explained by covariations and interplay between the rates of fluid-driven and collective entrainment of particles, as well as rates of particle deposition (Figs. 3-4) (see video supplement-Chartrand (2017)). At time scales much larger than bedload particle travel times of fractions or a few seconds or more (cf. Drake et al., 1988;Roseberry et al., 2012;Fathel et al., 2015), we envision that fluid-driven particle entrainment of any grain size fraction λ j sets a baseline 25 local activity, which can preferentially amplify the lower frequency content for smaller grain size classes (Figs. 3,4,12 and 13). The baseline flux is represented by a combination of a Poisson process for the transport behaviour of larger grain size classes, whereas smaller size classes exhibit a more time-persistent transport behaviour, with increasing correlation or memory at short time scales (Figs. 14 and 15) (Saletti et al., 2015).
Collective entrainment can further amplify the baseline activity across all overlapping frequencies, with a clear emphasis 30 on lower frequencies (Figs. 3 and 4). An emphasis on lower frequencies is related to the low-pass behaviour of Eqs. (8 and [15][16][17]. The amplifying effect of collective entrainment in setting the magnitude of the power spectra is maximized as µ j → σ j . Conversely, the collective entrainment effect diminishes for σ j >> µ j . Although we cannot definitely demonstrate that collective entrainment contributed to the grain size specific power spectra for PRE1 (Figs. 12 and 13)), two additional pieces of evidence support this interpretation. First, the PRE1 power spectra structure is generally consistent with that for our exploratory 3-particle birth-death system. Consistency occurs when the collective entrainment rate across all particles sizes is relatively large, and/or the autoregressive temporal correlation for smaller particle sizes is large relative to coarser particles (Fig. 4).
Second, and as mentioned above, histograms of particle flux for smaller particles sizes of PRE1 show agreement with the negative binomial distribution pmf (Figs. 9 and 11). This outcome supports the theoretical prediction of Ancey et al. (2008) for conditions when µ > 0, and the combination of evidence suggests an important role for collective entrainment to the transport of gravel mixtures (Fig. 3) (e.g. Ancey et al., 2006Ancey et al., , 2008Singh et al., 2009;Ancey, 2010;Heyman et al., 2013;Masteller and Finnegan, 2017;Lee and Jerolmack, 2018). As a result, our work adds to the growing recognition for the need for an expanded model of bedload transport (e.g. Ancey, 2010;Furbish et al., 2012a;Lee and Jerolmack, 2018;Yager et al., 5 2018b).
The local mixing of fluid and collective entrainment signals, less the effect of deposition, involves the physical coincidence of differing distributions of particle travel lengths and durations, and with particles at different points along a travel pathway (e.g. Fig. 3). In particular, we hypothesize that it is the interplay between these different phenomena combined which yields power spectra that fall off with increasing frequency for finer grain sizes, and which reduces to a white noise for grain sizes that 10 approach the mixture limiting size (Figs. 3,4,12 and 13). These results are not unique, and appear to represent a basic attribute of granular transport systems. For example, experimental results reported by Jerolmack and Paola (2010) and Dhont and Ancey (2018) show total transport power spectra which fall off at slopes between 1/f − 3.5 3 -1/f − 4.5 3 for increasing frequency, respectively. In Fig. C1 we provide composite spectral density distributions for PRE1, and find that fall off slopes for both steady-state and transient conditions range from 1/f − 5 3 -1/f − 4 3 for the frequency range 10 −3 -10 −0.5 . The frequency range reported by 15 Dhont and Ancey (2018) is notably lower than that reported here. Nonetheless, overall consistency between the three data sets in terms of the power spectra structure implies that similar bedload transport processes are responsible. Consequently, our work re-emphasizes the need to better understand controls on differences in how relatively small versus large grains are transported in rivers, because the mechanics of sediment mixture transport is not deterministic at any rate of motion (e.g. Einstein, 1937Einstein, , 1950Jackson and Beschta, 1982;Hassan et al., 1991;Ancey, 2010;Lajeunesse et al., 2010;Furbish et al., 2012a, b;Heyman 20 et al., 2013;Fathel et al., 2015;Saletti et al., 2015;Furbish et al., 2017;Masteller and Finnegan, 2017;Lee and Jerolmack, 2018;Yager et al., 2018a).

Concluding Remarks and Next Steps
Four time series of particle flux for steady-state and transient bed elevation reveals both consistent and contrasting transport behaviours across six different grain sizes. First, particle transport as bedload is characterized by a fluctuating time signal, 25 regardless of grain size. Second, the fluctuating signals encode transport behaviour differences that depend on grain size.
Notably, the motions of finer grain sizes are strongly correlated in time, whereas the transport of coarser grain sizes is weakly correlated in time. These results are supported by descriptive statistics, power spectral analysis and autocorrelation results, and are consistent with the results of Saletti et al. (2015).
An expanded Markov birth-death model (Ancey et al., 2008;Ancey, 2010;) offers a plausible explanation for 30 grain size dependent transport behaviour of gravel mixtures. The interplay between fluid-driven and collective entrainment, less the effect of deposition, locally amplifies the activity of increasingly smaller grain sizes. This affect is propagated downstream through the full distribution of particle motions, and the result is to increase the amplitude of the lower frequency component Autocorrelation Coefficient TIme Lag (seconds) Figure 15. Autocorrelation of windowed and detrended sediment particle flux time series for t2150 (gray curves, solid lines) and t4310 (blue curves, dashed lines).
of flux signals, with an emphasis on smaller grain sizes. Grain size dependent negative binomial fits to flux histograms further supports our interpretation, as well as the theoretical result of Ancey et al. (2008), and suggests collective entrainment was an important transport phenomenon in our experiments for the smaller grain sizes, likely due to bed impacts by larger grains.
Poisson distribution fits to flux histograms for the largest grain size classes suggests collective entrainment was not an important transport mechanism. In a more general way, we show that this overall cascading effect from large to smaller grains can be initiated in the birth-death model with an initial input signal of white noise, which could represent the motions of relatively large particles, or a system close to the transport threshold. Despite similarity between our experimental results and the expanded birth-death model, the link made here is conceptual, and substantial work remains.
For example, high resolution physical or possibly numerical experiments progressing from uniform to mixtures of several 5 size classes could be used to directly test and examine links between entrainment mechanisms and fluctuations of particle flux. Small flume setups would likely be ideal, and will build directly on previous work (Ancey et al., 2008;Ancey, 2010;Heyman et al., 2014;Lee and Jerolmack, 2018) by better identifying the parameter space for the fluid entrainment rate λ j , the rate constant of particle deposition σ j and the rate constant of collective entrainment for sediment mixtures µ jj . Furthermore, particle deposition is a key factor in the transport problem and the birth-death formulation, because it contributes directly to 10 the gain of a low pass filter which modulates the flux signal of differing grain sizes. However, our current understanding of particle deposition mechanics is insufficient to precisely quantify this effect in the experiments.
Code and data availability. The data and Jupyter Notebooks constructed to explore the three particle size birth-death model developed and proposed here, as well as analyze and compile the experimental results can be freely accessed at: https://doi.org/10.6084/m9.figshare.12268727 (Chartrand et al., 2020).
Video supplement. We provide two supplemental videos for download which illustrate sediment transport during a preliminary test experiment (Chartrand, 2017). The video options are at 6 Hz and 150 Hz, compared to the raw video collected at 30 Hz. The videos are useful for two reasons. First, the videos illustrate the general intensity and character of sediment transport discussed and reported here. Second, a collective entrainment event is captured and illustrated by both videos. The collective entrainment event represents a sequence of particle position adjustments which ultimately leads to a short term spike in local particle activity. Beginning at time ∼3 s of the 6 Hz video, the position of a large green particle in the frame center adjusts either due to fluid driven momentum delivery, or a particle collision that is difficult to make out due to depth of field limitations. This initial particle position adjustment leads to a small short term increase in local particle activity. At time ∼25 s, position adjustment of the same large green particle begins again, and about 15 s later a collective entrainment event begins, and there is a larger local short term increase of the local particle activity. The collective entrainment events are interesting because it is not clear if they are due to fluid driven phenomenon, or particle collisions. It is difficult to tell.

25
Appendix A: Average particle count The average particle count for each time series was calculated for different aggregation intervals (Figs. A1 and A2). In general, particle sizes smaller than 8 mm exhibit averages that vary among the different aggregation intervals, tending to steady values by the end of each time series. Grain sizes coarser than 8 mm, on the other hand, exhibit approximate steadiness across all aggregation intervals. Differences between the grain sizes in terms of the temporal trends of the aggregated average particle

Appendix B: Power Spectral Density for mass-based flux
The power spectral density for mass-based flux of PRE1 was calculated for all grain sizes and time series following the steps discussed above in Section 3.3 (Figs. B1 and B2). The mass-based power spectral density for each time series indicates that the 11.3 or 16 mm grain size classes have the largest spectral power. Overall, spectra structure is similar to that discussed for the count-based power spectra (Figs. 12 and 13). The finer grain size classes exhibit spectra that decrease exponentially with 5 increasing frequency proportional to 1/f − 6 3 to 1/f − 4 3 . The spectra fall off slopes diminish and approach white noise as grain sizes increase to the limiting maximum grain size. The primary difference between the mass-based and count-based spectral densities is the distribution of power among the grain size classes.
Appendix C: Composite Power Spectral Density for particle-based flux The power spectral density for the particle-based flux summed across all six grain size classes of PRE1 was calculated for the 10 transient time series following the steps discussed above in Section 3.3 (Fig. C1). The composite spectral density structure exhibits power that decreases exponentially with increasing frequency proportional to 1/f − 6 3 to 1/f − 4 3 , similar to that discussed for the count-based power spectra finer grain size classes (Figs. 12 and 13). Log normalized spectral density [-] PRE1t2150 and t4310min t1910min t4070min t2150min t4310min Figure C1. Estimated power spectral density for the total particle flux across all grain sizes for both the steady-state (t1910 and t4070) and transient (t2150 and t4310) time series. Each spectral density distribution has been normalized by the maximum power of the distribution.
The dashed red lines provide reference lines of different slope.

Appendix D: Power Spectra of the number of moving particles
The number of active particles n j of size class j for a specified streambed area A may be expressed as (Ancey, 2010): The variable v j is the immigration rate, λ j is the individual entrainment rate commonly attributed to fluid forces, µ j is a rate constant of collective entrainment, σ j is a rate constant of deposition, and w j is a rate constant of emigration. The Fourier 5 transform for the left hand side of Eq. D1 is found by differentiating the inverse transform with respect to t. Begin with the inverse Fourier transform: Differentiate both sides with respect to t: 10 The derivative on the right hand side is moved inside the integral because it operates on the exponential: Therefore, and, we recover the well known result We apply this result to the left hand side of Eq. D1, let a circumflex denote a Fourier transform, and ignore the (ω) dependence for brevity: where ω = 2π/T is the angular frequency with period T , and we have assumed thatν j ≈ŵ jnj for adjacent areas of the streambed A under statistically uniform transport. Rearranging, We specify A = σ − µ, and we write Eq. D8 in terms of its real and imaginary parts: The magnitude of the active number of particles in motion over A is: |n j | = A 2 + ω 2 (A 2 + ω 2 ) 2 0.5 |λ j |, |n j | = 1 (A 2 + ω 2 ) 0.5 |λ j |. (D10c) Using the definition of A, and noting that |n j | and |λ j | are just the signal amplitudes: Author contributions. The work reported here was an intellectual co-conspiracy between SMC and DJF-development of the three particle size birth-death model, the time series and Fourier analyses, preparation of the manuscript, etc. were all done together.
Competing interests. The authors declare no competing interests in preparation of this manuscript.