Dynamics and mechanics of bed-load tracer particles

Understanding the mechanics of bed load at the flood scale is necessary to link hydrology to landscape evolution. Here we report on observations of the transport of coarse sediment tracer particles in a cobblebedded alluvial river and a step-pool bedrock tributary, at the individual flood and multi-annual timescales. Tracer particle data for each survey are composed of measured displacement lengths for individual particles, and the number of tagged particles mobilized. For single floods we find that measured tracer particle displacement lengths are exponentially distributed; the number of mobile particles increases linearly with peak flood Shields stress, indicating partial bed load transport for all observed floods; and modal displacement distances scale linearly with excess shear velocity. These findings provide quantitative field support for a recently proposed modeling framework based on momentum conservation at the grain scale. Tracer displacement is weakly negatively correlated with particle size at the individual flood scale; however cumulative travel distance begins to show a stronger inverse relation to grain size when measured over many transport events. The observed spatial sorting of tracers approaches that of the river bed, and is consistent with size-selective deposition models and laboratory experiments. Tracer displacement data for the bedrock and alluvial channels collapse onto a single curve – despite more than an order of magnitude difference in channel slope – when variations of critical Shields stress and flow resistance between the two are accounted for. Results show how bed load dynamics may be predicted from a record of river stage, providing a direct link between climate and sediment transport.


Introduction
Understanding landscape denudation and its relation to climate requires an understanding of how a flood hydrograph drives the sediment mass flux leaving the system through rivers. While suspended sediment represents the largest fraction of mass exiting the landscape (Milliman and Syvitski, 1992;Willenbring et al., 2013), it is coarse bed load transport that sets the limiting rate of landscape incision through its control on bedrock erosion and channel geometry in gravel rivers (Sklar and Dietrich, 2004;Snyder et al., 2003;Parker et al., 2007). The rate of bed load transport is known to vary both spatially and temporally due to turbulence and granular phenomena such as clustering, bed forms, bed compaction, grain protrusion/hiding, and collective motion (Gomez, 1991;Kirchner et al., 1990;Schmeeckle et al., 2001;Strom et al., 2004;Ancey et al., 2008;Zimmermann et al., 2010;Marquis and Roy, 2012;Heyman et al., 2013), which makes predictions difficult (Recking et al., 2012) and point measurements highly variable (e.g., Gray et al., 2010). Bed load is especially difficult to predict near the threshold of motion (Recking et al., 2012), where transport is highly intermittent, often resulting in partial bed load transport, in which only a fraction of the bed is mobilized during a transporting event (Wilcock and McArdell, 1997). Further confounding predictions is that many gravel streams adjust their geometry to an effective discharge (Wolman and Miller, 1960), which occurs at a flow slightly above (1.2-1.4 times) the threshold of motion for the median grain size (Parker, 1978;Parker et al., , 2007, indicating that partial transport may be the dominant transport regime within gravel rivers. The spatially variable and highly intermittent flux during partial transport (Wilcock and McArdell, 1997;Haschenburger and Wilcock, 2003), compounded with the added difficulty of a varying Published by Copernicus Publications on behalf of the European Geosciences Union. 514 C. B. Phillips and D. J. Jerolmack: Dynamics and mechanics of bed-load tracer particles sediment supply, necessitates long-term observations to decipher bed load dynamics in the field.
Passive tracer particles, in particular passive integrated transponder radio-frequency identification (PIT RFID) tagged particles, are becoming an attractive low-cost and low-maintenance method of measuring bed load particle dynamics. The application of passive tracer particles has taken various forms, such as exotic lithologies (Houbrechts et al., 2011), painted bed material (Wilcock, 1997b), magnetic (Hassan et al., 1991), radioactive (Sayre and Hubbell, 1965;Bradley et al., 2010), and RFID (Lamarre et al., 2005;Bradley and Tucker, 2012;Phillips et al., 2013;Schneider et al., 2014). A benefit of RFID-equipped tracer particles is that each particle is uniquely identified, which allows its position to be measured at longer timescales with high recovery rates (Bradley and Tucker, 2012;Phillips et al., 2013). An advantage of long-term (multi-flood to multi-annual) observations of tracer particles is that they sample over temporal variations in fluid stress and spatial heterogeneity in the river bed, and thus present an integrated picture of bed load transport dynamics.
In this paper we present the results of a 2-year deployment of several populations of RFID tracer cobbles within alluvial and bedrock sections of a river, for single-flood and yearly timescales. Throughout this manuscript, we define bedrock channels following the definition of Turowski et al. (2008) as channels that cannot substantially alter their geometry without eroding bedrock; conversely, alluvial channels are channels that are able to freely adjust their geometry without eroding bedrock. At the individual flood scale we examine the tracer displacement distributions, and the fraction of tracers mobilized. We show that tracer displacements and the fraction mobile are consistent with results from a recent momentum conservation framework (Charru et al., 2004;Lajeunesse et al., 2010). We employ a recently developed dimensionless impulse framework (Phillips et al., 2013) to account for unsteadiness in the hydrograph, which allows us to apply a fluid momentum conservation approach to long-term tracer displacement data. For flows within the partial transport regime, we demonstrate that tracer displacements are short and close to the limit of one step per flood. Furthermore, we demonstrate the generality of the long-term tracer displacement results within the main channel with a smaller deployment of tracers in a step-pool bedrock tributary. We show that, by accounting for flow resistance and differences in the threshold of motion, displacement dynamics in the step-pool and main channels are similar. Lastly, we analyze and compare the sorting of tracer particles with that of the river to show that the emerging sorting patterns are consistent with a sizeselective transport sorting model.

Theory
In the following sections we present the relevant theoretical background that guides the analysis and interpretations of our tracer particle results. The theoretical background is intended as a brief introduction to the topics of sediment transport mechanics and dynamics, quantifying hydrologic forcing, and the downstream sorting of sediment by particle size.

Sediment transport at the particle scale
Under a wide range of bed load transport conditions, coarse sediment particles undergo short steps separated by longer periods of rest, which leads to probabilistic descriptions of particle motion (Einstein, 1937). A particle step is defined as the distance the particle is transported from entrainment to deposition, and the rest duration is the time between deposition and subsequent entrainment. The motion of particles in bed load transport is comprised of sliding, rolling, or short hops called saltations (Drake et al., 1988), where the travel time is generally much shorter than the rest duration (Lajeunesse et al., 2010;Martin et al., 2012;Furbish et al., 2012b, a;Roseberry et al., 2012). For near-threshold bed load transport, in which only bed surface particles are mobile, bed load flux may be described as the product of the particle velocity and surface density (particles/area) of moving grains (Bridge and Dominic, 1984;Wiberg and Smith, 1989;Parker et al., 2003;Lajeunesse et al., 2010;Furbish et al., 2012b), or similarly the product of the particle entrainment rate and the average particle step length (Einstein, 1950;Wilcock, 1997a;Wong et al., 2007;Ganti et al., 2010;Furbish et al., 2012b). The combination of particle velocity, number of mobile surface particles, depth of the mobile layer, and a threshold stress typically result in a nonlinear relationship between bed load flux and the fluid shear stress (Meyer-Petter and Muller, 1948;Fernandez Luque and Van Beek, 1976;Wong and Parker, 2006;Furbish et al., 2012b). Here it should be noted that the above formulations for bed load particle flux both require averaging the measured quantities over yet-undetermined timescales (Ancey, 2010;Furbish et al., 2012b). For steady turbulent flows in the laboratory, the particle velocity and step length have been shown to scale linearly with the excess shear velocity (U * − U * c ) (Fernandez Luque and Van Beek, 1976;Lajeunesse et al., 2010;Roseberry et al., 2012;Martin et al., 2012), where U * = √ ρτ b is the shear velocity (m s −1 ), τ b is the basal shear stress, ρ is the fluid density (1000 kg m −3 ), and U * c (m s −1 ) is the threshold shear velocity for initiation of sediment motion. Specifically, Lajeunesse et al. (2010) found that the modal particle step length scales as where X is the transport distance (m), D is the particle median axis (m), C = 70 is an empirically determined constant, and V s = √ R g D is the settling velocity in the limit of large Earth Surf. Dynam., 2, 513-530, 2014 www.earth-surf-dynam.net/2/513/2014/ particle Reynolds numbers, where R is the submerged specific gravity of the particles and g is the acceleration due to gravity (m s −2 ). The surface density of moving particles was found to increase linearly with the Shields stress (τ * ) (Lajeunesse et al., 2010), where τ * = τ b /(ρ s − ρ)g D, and ρ s is the sediment density (2650 kg m −3 for quartz). That the dependencies of the step length, particle velocity, and mobile surface density on shear velocity have been recently validated for both unimodal and bimodal grain size distributions under turbulent flow (Lajeunesse et al., 2010;Houssais and Lajeunesse, 2012) encourages us to extend these results to interpret tracer particle data at the field scale.
Treating the particle behavior probabilistically, we focus on the distributions of particle steps and rests. In the laboratory, the distribution of particle step lengths for a given stress and grain size have been observed to follow exponential or gamma-like distributions (Lajeunesse et al., 2010;Hill et al., 2010;Martin et al., 2012;Roseberry et al., 2012). However, for mixed grain size distributions, heavy-tailed statistics can emerge due to a summation of exponential step lengths for each size group (Hill et al., 2010). Examining passive tracers in the field introduces an ambiguity; one measures particle displacement -i.e., the distance a particle travels between successive surveys of its position -but this displacement is composed of an unknown number of steps and rests. Displacement length distributions measured for individual floods, and at longer timescales over many floods, typically follow exponential or gamma-like distributions (Hassan et al., 1991;Schmidt and Ergenzinger, 1992;Habersack, 2001;Lamarre and Roy, 2008;Bradley and Tucker, 2012;Hassan et al., 2013;Phillips et al., 2013). We propose two simple limits for particle displacement during a flood: (1) the lower limit is that a particle executes a single step, with a characteristic length scale predicted by Eq. (1), and (2) the upper limit is continuous particle transport, with no rests, for the duration of the flow that exceeds the threshold entrainment stress. We explore tracer displacements within the context of these two limits.
Upon deposition, the rest duration before subsequent entrainment is constrained by two criteria: first the stress must exceed the threshold of motion locally, and second the particle must be exposed to the flow (Martin et al., 2012). The stochastic erosion and deposition of the river-bed surface acts to bury and excavate particles, and recent laboratory results suggest that this produces heavy-tailed particle rest durations (Martin et al., 2012(Martin et al., , 2014. Although these rest durations cannot be measured from passive tracers in the field, our previous work used the dispersion of the tracer plume to infer similar behavior to laboratory experiments (Phillips et al., 2013). Accordingly, we will not consider the particle rest duration or tracer dispersion in this article.

Dimensionless impulse
Flows in natural coarse-grained rivers are inherently unsteady: from the microscopic scale of variations in turbulence, to macroscopic fluctuations in discharge within a flood, to the rise and fall of the hydrograph throughout a series of floods. At the smallest relevant scales of turbulence, the threshold of motion is determined by the impulse, the product of shear stress magnitude, and duration (Diplas et al., 2008). Due in part to the difficulties in measuring tracer particle motion and near-bed stresses during floods, the fluid shear stress is commonly quantified through use of a bulkflow parameter such as the depth-slope product, τ b = ρg h S (Church and Hassan, 1992;Hassan et al., 1991;Ferguson and Wathen, 1998;Haschenburger and Church, 1998;Lenzi, 2004;Haschenburger, 2011), where h is the flow depth (m) and S is channel slope. For coarse-grained streams this simplification is perhaps more reasonable, as particle inertial timescales are large and thus coarse particles are insensitive to a range of turbulent stress fluctuations (Diplas et al., 2008;Celik et al., 2010;Valyrakis et al., 2010Valyrakis et al., , 2013. Although some readers may object to the assumption of steady and uniform flow for a flood, the flow may be considered quasisteady so long as the hydrograph varies slowly compared to the grain inertial timescale (on the order of several seconds), and quasi-uniform so long as water surface slope remains approximately constant. Thus for large ensembles of particles over many floods, the idea of employing a normal flow approximation becomes tenable. Accordingly, Phillips et al. (2013) introduced the dimensionless impulse to quantify the time-integrated fluid momentum in excess of threshold, assuming normal flow. Here D 50 represents the median grain size of the tracers, t s represents the start of a flood, and t f represents the end of a flood of interest. The integral is only calculated over the record of U * > U * c , as subthreshold flows do not transport sediment. We note that I * in this study represents a cumulative metric of reach-averaged fluid momentum.

Downstream sediment sorting
The spatial pattern of diminishing grain size with increasing distance from the headwaters is near universal among gravel rivers, and results from a combination of size-selective sorting and particle abrasion (e.g., Paola et al., 1992;Kodama, 1994;Paola and Seal, 1995;Ferguson et al., 1996;Gasparini et al., 1999Gasparini et al., , 2004Fedele and Paola, 2007;Jerolmack and Brzinski, 2010). For tracer particles, the relatively short distances and timescales involved preclude abrasion as a mechanism for observed downstream fining (Ferguson et al., 1996). Thus we further explore the mechanisms and implications 516 C. B. Phillips and D. J. Jerolmack: Dynamics and mechanics of bed-load tracer particles of size-selective transport as it pertains to tracer particles and the river bed. We look to laboratory experiments to inform the following analysis, as it is difficult to generalize the rate at which tracer particles sort from previous field studies. Flume experiments with a heterogeneous sediment input show that particles rapidly segregate by size to achieve an equilibrium profile, and that subsequent transport results in an elongation (stretching) of this sorting profile (Paola et al., 1992;Paola and Seal, 1995;Seal et al., 1997;Toro-Escobar et al., 2000). The self-similar sorting profile means that longitudinal profiles collapse onto a single curve when downstream distance (X) is cast as a dimensionless extraction length X * = X/L, where L (m) is the distance at which 100 % of the coarse material in transport has been extracted (deposited). In laboratory experiments and in natural rivers, L is taken as the distance from the input/source to the gravelsand transition (Paola and Seal, 1995); X * = 1 indicates that all gravel particles are deposited upstream of this location.
In the absence of a well-defined gravel front the extraction length can be difficult to determine, though Toro-Escobar et al. (2000) suggest that the location where 90 or 95 % of the source material has been extracted (deposited) is a suitable proxy for L. In the case of tracer particles, this represents the distance downstream from the source to the location of 95 % recovery. Laboratory and modeling results (Fedele and Paola, 2007) demonstrate that the downstream sorting of a gravel mixture can be described by its initial variance and mean in the following formulations: where σ o and D o are the standard deviation and mean of the input material, D is the mean grain size, σ is the standard deviation of the grain size distribution, and C 1 and C 2 are constants. For selective deposition in gravel rivers, D and σ should decrease exponentially at approximately the same proportion, resulting in a constant coefficient of variation (σ/D) (Fedele and Paola, 2007). These results also suggest that, due to the self-similar sorting profile, separate populations of tracer particles with similar initial grain size distributions (e.g., D and σ ) should exhibit similar dynamics when properly rescaled by X * .

Field site
Field deployment of tracer particles took place in the Mameyes River basin, located within the Luquillo Critical Zone Observatory in northeastern Puerto Rico. Coarsegrained tracers equipped with PIT RFID tags were deployed in the main stem of the Mameyes River, and in a steep tributary (Fig. 1a). The Mameyes River drains the center of the Luquillo Mountain Range, which commonly receives more than 4000 mm yr −1 of precipitation in the headwaters due to a large orographic effect. Precipitation occurs frequently throughout the year in the form of high-intensity, short-duration events, resulting in frequent flash flooding (Schellekens et al., 2004). The 1.2 km study reach in the Mameyes begins just downstream of where the river exits the mountains (drainage area of 24.21 km 2 ). At the main channel field site, stage was recorded at 5 min intervals for 40 consecutive days by an In-Situ Level Troll 500 and measured from surveys of high flow debris following the largest floods, which were correlated to discharge (Q) measured 3.5 km upstream by a US Geological Survey (USGS) gage (gage no. 50065500, 15 min resolution) to obtain a reach-averaged depth record for the study period (Fig. 2a). The combination of the automated stage measurements and high flow surveys captures the entire range of discharge for the 2-year study period. For calculations of the frictional resistance we use the hydraulic radius (h r ). This section of the Mameyes River ( Fig. 3a) exhibits minimal meandering with nearly constant width (20 m), and has a slope of S = 7.8 × 10 −3 . The slope was extracted from a lidar digital elevation model (DEM) along the channel center (1 m horizontal and vertical resolution) (Fig. 3b). The smaller headwaters tributary is located in catchment three (drainage area of 0.58 km 2 ) of the Bisley Experimental Watershed (Bisley 3). The Bisley 3 study reach is characterized as a step-pool stream with a slope of S = 1.2 × 10 −1 , width ranging from 2 to 4 m, and boulder steps that range from 0.5 to 2 m in height (Fig. 3c). For this reach ( Fig. 3c) the slope was determined from a longitudinal profile surveyed in the field (Fig. 3d); due to extremely dense forest canopy, a longitudinal profile extracted from the lidar DEM does not accurately represent the heterogeneity in channel topography (Fig. 3c). It should be noted that a simple linear regression is unlikely to capture the heterogeneity in transport slopes for this step-pool stream; however we use it to remain consistent with our main channel field site. Stage was recorded at 1 min intervals at the field site for 59 consecutive days by an In-Situ Level Troll 500, which was correlated to discharge measured by a US Forest Service (USFS) gage (15 min resolution) located ∼ 100 m upstream to obtain a reach-averaged depth record for the duration of the study (Fig. 2b).
To characterize the spatial sorting of the stream bed, we performed pebble counts at 200 m intervals from the start of the Mameyes RFID tracers to the perceived gravel-sand transition. When analyzing the stream sorting we restrict the analysis to the depositional part of the river, and thus we only use measurements downstream of the start of the alluvial plain (approximately 200 m downstream of tracer installation location). The D and σ of the tracers were measured at the center of eight linearly spaced bins moving downstream from the initial placement location. The number of bins was determined to ensure enough tracers within each bin for accurate statistics. To determine the extraction length X * for the stream, we set the basin length (L) as the distance from the start of the alluvial plain to the perceived gravel-sand transition. The exact location of the gravel-sand transition in the Mameyes River is obscured in the field due to substantial anthropogenic modification of the river. However, with the use of airborne lidar data we have identified a substantial break in channel slope within the confines of the golf course approximately 3.5 km downstream of the start of the alluvial plain, which we take as the gravel-sand transition. As there is no well-defined front for both populations of tracer particles, we set L as the distance at which 95 % of the tracers recovered remain upstream for each population of tracers at the time of the final survey. Ideally L is the distance at which 100 % of the total (recovered and unrecoverable) tracers remain upstream; however we cannot accurately determine this point. In order to reduce the variability in the surface grain size distributions we determined the D 84 from an exponential fit to pebble count data collected up-and downstream of each study reach.

Tracer particles
Two populations of 150 tracers were installed in the summers of 2010 and 2011 in the Mameyes as it exits the mountains near the start of the alluvial plain. A smaller population of 51 tracers was installed in the Bisley 3 stream in the summer of 2010. All three tracer particle populations are composed of cobbles from the stream bed and have narrow grain size distributions centered on the bed D 50 (Fig. 1c and d). Narrow grain size distributions were selected to facilitate equal mobility (e.g., Wiberg and Smith, 1987) within the tracer populations. The median grain size values for both tracer populations and the river bed at the Mameyes site were 12, 13, and 11 cm, respectively. Median particle diameters for the Bisley 3 stream and tracer particles were 12 and 13.5 cm, respectively. The median grain size for the main channel reach represents the average of three separate pebble counts of 100 particles each (Wolman, 1954). In the Mameyes and Bisley 3 field sites the tracers are fully submerged (average at both sites h/D 50 = 7) during transport for flows above the threshold of motion. Both populations of tracers in the Mameyes River were deployed in the same reach (Fig. 1e) in a 20 m × 20 m grid with 1 m spacing across the width of  (1) were 100, 99, and 94.6 %. The initial recovery rate for population 1 is low as a result of an incomplete survey, which was cut short for safety concerns due to the occurrence of a second large flood. The second population of tracers was placed in the river as two installments of 80 and 70 tracers on two consecutive days due to a small flood, which resulted in a minor amount of burial from fine sediment for the initial 80 tracers installed. Impacts of the initial increased embeddedness on this subset of tracers are not observable at the multi-flood scale. Unrecovered tracers have the potential to bias the mean value, as it is possible that tracers could be buried beyond the detection limit, destroyed, missed, or are farther downstream. However, it was   Fig. 1b and f. Surveyed positions for both field sites were transformed from Cartesian coordinates to a stream-wise normal system using a methodology similar to that developed by Legleiter and Kyriakidis (2007). Tracer particles were located using two wands (manufactured by Oregon RFID) with empirically determined horizontal detection limits of 50 and 20 cm, respectively. PIT RFID tags (32 mm HDX tags supplied by Oregon RFID) were detectable when buried at depths up to 10-20 cm below the river bed for the small wand (depending on tag orientation), and 50 cm for the large wand. The maximum combined survey and detection error is estimated to be 1 m and 45 cm for the large and small wands, respectively (for a discussion of wand detection distances and limitations see Chapuis et al., 2014). All measured tracer motion recorded below the detection threshold was considered to be error and set to zero. Assuming that tracers detectable with the larger wand only were buried beyond the detection limit of the small wand, then only 6 % of the tracers were buried more than 20 cm below the surface. This suggests that the majority of the tracers were near the surface. All tracer data are available online as part of the Luquillo Critical Zone Observatory database (see acknowledgements).

Hydrologic forcing
For the study reaches, U * = √ g h S was estimated assuming steady and uniform flow; Shields stress was also estimated for comparison to other studies. Long-term flow records for the Bisley 3 stream are measured near a series of large culverts that artificially truncate the largest floods; this does not affect the calculations of I * greatly, but does add ambiguity to the distributions examined in following sections. Therefore, we limit our analysis to the hydrograph on the main channel of the Mameyes River, as it represents the highest quality data.
The value of I * (Eq. 2) was computed for each flood above the threshold of motion (Fig. 2), and also for the cumulative time periods between successive tracer surveys. The calculation of I * is particularly sensitive to the value of U * c , a parameter that is known to vary both temporally and spatially within a flood and from flood to flood (Kirchner et al., 1990;Charru et al., 2004;Turowski et al., 2011;Marquis and Roy, 2012). We treat U * c as a constant by necessity as we lack a theoretical or empirical methodology with which to account for these effects. We determined the value of U * c in two manners: (1) U * c1 was calculated from the fraction of mobile tracers for individual floods, and (2) U * c2 was determined as the value that provided the best collapse of the mean tracer displacement data (both methods are used in Sects. 4.2.1 and 4.2.2). Using two definitions of critical shear velocity, the resulting values of I * for a representative flood are 3815 for U * c1 and 806 for U * c2 (Fig. 4). The discrepancy between the values of I * is due to the broadening asymmetric shape of the hydrograph, where a small change in U * c  Fig. 2a. The limits of integration for U * c1 and U * c2 are t 1 to t 4 and t 2 to t 3 , respectively. The shaded region represents the region integrated for the calculation of I * using U * c1 .
can result in an order of magnitude increase in I * . When calculated over numerous floods, uncertainty in U * c can result in substantial differences in I * . However, it should be noted that both approaches produce the same scaling relationships, but differ in the magnitude of the coefficients. A potential drawback of calculating U * and U * c from a reach-averaged depth-slope product is that large instantaneous values of U * due to turbulent fluctuations cannot be accounted for. This simplification could result in tracer movement from turbulent fluctuations being attributed to a reach-averaged stress, and hence a biased estimate for the actual threshold value. The upshot is that estimates for I * are least accurate for lowmagnitude, short-duration floods. Data points likely to be affected by this consequence of determining a reach-averaged U * c are included for completeness and indicated in the following figures where appropriate. Due to this drawback we do not recommend using Eq.
(2) for individual floods without an independent measure of U * c . For the remainder of the manuscript, except where noted, we use the value of U * c as determined by the fraction of mobile tracers (U * c1 ).

Hydrology filtered through sediment mechanics and dimensionless impulse
River discharge is the most commonly reported variable in relating long-term sediment dynamics to hydrologic forcing. However, the momentum framework presented in Sect. 2.1 reminds us that fluid stress -rather than water discharge -is the relevant parameter to consider for driving sediment motion. For considerations of bed load transport, the threshold of motion applies a filter to these data; only flows that exceed the critical stress for entrainment are relevant for assessing particle transport. Accordingly, we empirically estimate the threshold of motion by determining the fraction of mobi- lized tracers (f ) for several individual floods, where tracers' positions were surveyed before and after the event. Based on a momentum balance approach (cf. Lajeunesse et al., 2010), we anticipate that f should scale linearly with Shields stress (Fig. 5). A linear relation provides a reasonable fit to the data when plotted against the peak Shields stress of the flood for the D 50 of the tracer population. We treat the intercept as the threshold of motion, determining that τ * c = 0.023 (U * c1 = 0.22 m s −1 ) (Fig. 5). We use the Shields stress at the flood peak -rather than average Shields stress -as it does not require additional information; computing the average stress associated with a flood requires choosing a threshold value for fluid stress. Error bars for f account for the number of missing tracers, where the upper and lower lines indicate the absolute maximum and minimum values for the fraction mobile by assuming that all tracers not recovered moved or did not move, respectively. For all flood events monitored, the fraction mobile remained well below 1 (f < 1). The magnitude-frequency distribution for U * for the entire period of record (Fig. 6a) presents a fuller picture of the statistical scaling of flow within the hydrograph. For small to intermediate values of U * the curve appears to be a straight line on a log-log plot (Fig. 6a), indicating potential power-law scaling for this region, a common feature of flood hydrology (Turcotte, 1994;Lague et al., 2005;Molnar et al., 2006). The power-law scaling is even more evident for discharge (Fig. 6a, inset). The magnitude-frequency distribution for both U * and Q exhibits a truncation to the power-law scaling that occurs at approximately the threshold of motion (Fig. 6a). The upper truncation of the distribution of U * is well fit by an exponential function, indicating that the shear velocity for flows exceeding threshold is well described by a single average value, U * = 0.27 m s −1 ( τ * = 0.033, U * /U * c1 = 1.23). The exponential decay in probability for shear velocities above critical further indicates that most bed load transport events in the Mameyes , where U * is the average value of the distribution of U * > U * c1 , t is the average duration of a flood above the threshold of motion for the duration of the study, and D 50 is the median grain size of the tracer particles.
do not exceed the threshold of motion by much. Despite a large range in discharge (Q peak /Q c = 14.6), the peak shear velocity U * peak = 0.46 m s −1 (τ * peak = 0.1) observed only exceeded U * c1 by a factor of 2. This filtering of the discharge data is particularly strong as the highest discharge values in the 20-year instantaneous record exceed our peak observed flood by a factor of 3.1 (Q = 444.6 m 3 s −1 ), while the highest shear velocity values associated with these extreme events are only 1.28 times (U * = 0.59 m s −1 ) as high as the peak observed value in this study.
We analyze the magnitude-frequency distribution of I * (Fig. 6b) and find that the distribution of I * is composed of several scaling regions, though it does not appear to be heavy tailed. The smallest values of I * are artificially truncated by the resolution of the river stage measurements (15 min). The probability density function (PDF) of I * (Fig. 6b, inset) has a pronounced peak that coincides with the product of the average excess shear velocity and average flood dura-tion, ( U * − U * c ) t /D 50 = 2070. Here U * is the average value of the magnitude frequency distribution of U * > U * c (Fig. 6a), and t is the average duration of a flood above the threshold of motion for the period of record. The distributions of U * and I * calculated over the relatively short study duration (2 years) are the same as those calculated using a longer flow record (20 years available from the USGS) for the Mameyes River. Thus, data indicate that there is a welldefined "characteristic flood" for the Mameyes.

Individual flood scale
The distribution of particle displacements resulting from several floods was determined from surveys of both populations of tracers at the Mameyes site, in the summers of 2010 and 2011. We normalize each tracer's transport distance (X i ) by its median diameter (D i ). For individual flooding events above the threshold of motion, the majority of the cumulative distribution functions (CDF) of tracer particle displacement are well described by exponential functions (Fig. 7a), except for two tracer displacement CDFs which decay faster than exponential functions (green circles and cyan squares in Fig. 7a). Each tracer displacement CDF is normalized by its mean displacement ( X/D ) to facilitate plotting all single events in one graph (Fig. 7a). Typical travel distances for individual tracers were on the order of a few meters for each flood. We compare the dimensionless tracer distance for all single events against the peak shear velocity for that event (Fig. 7b), normalized by V s , and find that the modal tracer displacement is well described by a linear relation, as anticipated by the momentum framework presented earlier. A fit to the modal displacement distances provides another estimate of the threshold stress, U * c = 0.13 (τ * c = 0.016), from the intercept. This is likely a lower estimate of U * c , and is not the same value as determined previously for the fraction mobile data. Finally, we plot Eq. (1) over a contour density map of our field data (Fig. 7b), using k = 70 from Lajeunesse et al. (2010). Remarkably, the modal step lengths predicted from the laboratory-derived relation of Lajeunesse et al. (2010) (Eq. 1) run through the modal displacement distances measured from our field data in the Mameyes.

Multi-flood scale
At the multi-flood scale we analyze the long-term behavior of the tracer particles' displacement. We normalize the CDFs of cumulative travel distance by each survey's mean value, which results in a collapse of the data. This collapse suggests that the mean value is a reasonable descriptor of the dynamics of each tracer population. We note here that the CDFs in Fig. 8 are truncated at the lower end due to measurement accuracy, and at the upper end of the distribution due to unrecovered tracers (see Hassan et al., 2013, for a discussion of the effects of unrecovered tracers on the scaling of Contour density plot of X/D against the excess shear velocity normalized by the settling velocity for each tracer. The contour colors represent the density of tracers within that location. The value of U * is for the flood peak, while the value of U * c is treated as a fitting parameter. The black line represents the expected linear relationship between the dimensionless shear velocity and the modal tracer step length (Eq. 1). the CDF, while for a discussion of the functional distribution for the CDF see Bradley and Tucker, 2012). We therefore analyze the mean cumulative tracer travel distances ( X/D ) using the dimensionless impulse (I * ) Eq. (2). We find that the X/D scales linearly with I * for both populations of tracer particles (Fig. 9). Due to a limited number of repeat surveys we utilize all permutations of tracer surveys; in other words, we determine the X/D and I * for all possible sampling intervals. Using all permutations of tracer surveys does require the assumption that the sequence of floods does not exert substantial control on the mechanics of particle displacement. However, flood sequence and particle embeddedness may explain some of the scatter in the X/D displacement data (Fig. 9). Due to the close agreement between our field data and laboratory results (Fig. 7b) we use Eq. (1) to calcu- late the two limits of particle transport discussed in Sect. 2.1. Limit 1, in which entrained particles execute one step per flood, was calculated using Eq. (1) with the peak U * for each flood for the study duration above U * c1 . Limit 2, continuous motion, was calculated as the product of Eq. (1) and t/T s , where t is the duration of the flood above the threshold of motion, and T s = 10.6 √ D 50 /R g is the expected particle step duration from Lajeunesse et al. (2010). For limit 2 we use the average shear velocity over the duration of a flood in Eq. (1). When these limits are compared with the tracer particle data, we find that the tracer particles' mean displacement is significantly closer to limit 1 than limit 2 (Fig. 9), consistent with highly intermittent and partial bed load transport.
Due to time limitations in the field we were unable to survey more than three individual floods for the Bisley 3 stream, and thus are unable to determine the threshold of motion from the fraction of mobile tracers. In order to compare the Mameyes data with the Bisley 3 tracer data, values for U * c from both sites must be determined using the same methodology (Wilcock, 1988). Therefore we utilize method two, described in Sect. 3, to determine the values of U * c2 for this comparison. The values of U * c2 that provide the best collapse of the tracer data are 0.28 m s −1 (τ * c2 = 0.038) and 0.81 m s −1 (τ * c2 = 0.303) for the Mameyes and Bisley 3 sites, respectively. We find that the Bisley 3 data are also well characterized by a linear relation, with a slope that is 1 order of magnitude lower than that of the Mameyes site (Fig. 10). It is intriguing that the Bisley 3 displacement data form a well-defined linear trend, considering the limited number of tracers and the rough channel geometry ( Fig. 3c and d). We note here that the slope of the linear relation computed with U * c2 for the Mameyes site is an order of magnitude larger than that determined using U * c1 ; however, the linear form of the relationship is robust for a wide range of threshold values.

Tracer particle sorting in the Mameyes
At longer timescales, tracer particle sorting by size is readily apparent in the Mameyes, as observed in other studies (Hassan et al., 1991;Ferguson et al., 1996;Hodge et al., 2011). Sorting may be the result of (1) an inverse relation between particle step length and grain size, and/or (2) differences in entrainment frequency throughout a flood as a function of grain size. To examine differences in entrainment, we separate the tracer particles into two populations: particles that moved at least once for all single floods, and particles that remained immobile. All particle size distributions are well fit by lognormal functions. We aggregate all of the mobile and immobile particles into two combined distributions for statistical purposes (Fig. 11a). The size distribution of mobile particles is finer at the coarse end of the distribution compared to immobile particles (Fig. 11a), though there is a fair amount of overlap between the CDFs. We used a two-sample t test (equal variance, unequal sample size) on the natural-log-transformed distributions to determine that the difference in the mean values ( mobile = 119 mm, immobile = 137 mm) between the mobile and immobile populations is statistically different (t statistic = 4.02, degrees freedom = 345, p value < 0.001). At the single-flood scale there does not appear to be a significant dependence of displacement length on particle size. When we compare the displacements during a flood to the expected step length calculated from Eq. (1), it becomes apparent that the majority of tracers have displacements that are close to the expected step length, while a significant number of smaller (D i < D 50 ) particles have longer displacements (Fig. 11b). The cumulative effects of minor grain size sorting at the flood scale (Fig. 11) result in the rapid development of tracer sorting at the annual scale for the Mameyes River. In order to connect the dynamics of tracers to the downstream fining pattern of river bed, we first analyze the downstream grain size trend of the bed surface of the Mameyes River as determined from pebble counts. The starting location of the tracers coincides with where the Mameyes River exits the mountains (at approximately 0 km downstream in Fig. 12a). The start of the alluvial plain begins at approximately 0.5-0.75 km downstream in Fig. 12a. From this point downstream the bed D 84 fines by roughly a factor of 2, while there is only minimal decrease in the bed D 50 and D 16 (Fig. 12a). This suggests that downstream fining of the river occurs through deposition of the coarsest particles. Starting at X = 300 m (Fig. 12a), the standard deviation and mean grain size decrease exponentially with initial statistics of σ o = 185 mm and D = 185 mm ( Fig. 12b and e). The exponent of the fitted exponential function (C 1 in Eq. 3) is 1.2 for the stream. This results in an approximately constant coefficient of variation σ/D = 0.93 (Fig. 12c). Using the parameters determined in Fig. 12b and c, we apply the full model of Fedele and Paola (2007) (Eq. 4) to the spatial decrease in D and find that the model seems to underpredict the mean value (Fig. 12e). Potential reasons for the underfit are given in Sect. 5.2.
Turning to the tracers emplaced in the Mameyes, sorting of both populations by size is readily apparent (Fig. 12a). The spatial decrease in σ for both tracer populations (σ P1 and σ P2 ) are well described by exponential functions (Fig. 12b). The exponents of the fitted exponential functions (C 1 in Eq. 3) for population 1 and 2 are 0.83 and 0.73, respectively. The initial standard deviations of tracer populations 1 and 2 at placement are 48 and 41 mm, respectively. These values are very close to the coefficients for the fitted exponential function in Fig. 12a. The initial mean values for tracer populations 1 and 2 are 129 and 140 mm, respectively. In accordance with the model of Fedele and Paola (2007) we find that both populations of tracers have nearly constant coefficients of variation (σ/D) of 0.25 and 0.26, respectively (Fig. 12c). The coefficient of variation for the tracers is expected to be low due to the narrow grain size distribution. Furthermore, we find that the σ data can be reasonably collapsed by normalizing by σ o from the upstream end of the depositional system (Fig. 12d). Here we apply Eq. (4) using the parameters determined in Fig. 12b and c and find that it underpredicts the rate at which D decays. In applying Eq. (4) for the tracers we set D o = 173 mm to the value for the stream at the location where they were placed. This is because the model is for a bed with a continuous source where the mean at X * = 0 does not change, while for a finite population of tracers the coarser particles are deposited, and thus over time D at X * = 0 will coarsen. Potential reasons for the underfit are given in Sect. 5.2.

Sediment mechanics
At the single-flood scale, tracer particle displacements have been shown to be well described by exponential, gamma, and power-law distributions (Phillips et al., 2013;Habersack, 2001;Hassan et al., 2013). Here we find that the majority of observations are well described by exponential distributions, with two exceptions that decay faster than exponentially (Fig. 7a). The distributions that decay faster than exponentially are likely due to undersampling as a result of the small number of mobile tracer particles in these floods. The particle velocity distribution has been shown in laboratory experiments and theoretically to scale exponentially, except for small sample sizes for which the distribution decays faster than exponential (Furbish and Schmeeckle, 2013); this pattern is also likely true for the distribution of particle displacements. Therefore, we might expect that a larger population of tracers would converge to the exponential distribution. The remarkable agreement between the modal tracer displacement data at the single-flood scale for the Mameyes and the laboratory results of Lajeunesse et al. (2010) (Fig. 7b) strongly suggests that the most likely tracer displacement for these floods is a single step length. This further suggests that the observed exponential distributions for tracer displacement (Fig. 7a) represent the distributions of tracer step lengths. The observation of mostly single step displacements holds for all observed individual floods, despite a 3fold increase in Shields stress (τ * /τ * c1 = 3.04), demonstrating that partial and intermittent bed load transport occurred under all observed conditions. This is consistent with previous field observations by Mao and Surian (2010), who reported that partial transport occurred over a 3-fold increase in Shields stress. However, tracer displacement length distributions become increasingly skewed toward larger values with increasing flood strength ( Fig. 7a and b), indicating that increasing numbers of particles experience more frequent re-entrainment as shear stress increases. Nonetheless, the fraction of mobile tracers remained significantly below 1 for all observed floods (Fig. 5), reinforcing the inference that partial bed load transport is the dominant mode of transport here. Plots of mean tracer displacement length against cumulative impulse show that tracer motion is far from the continuous limit, and much closer to the lower limit of one step per flood (Fig. 9). The linear scaling observed (Fig. 9) for both populations of Mameyes tracers indicates that, to first order, the mean displacement is determined by the total momentum imparted to the stream bed. The agreement in the two populations suggests that effects of embeddedness and initial placement are not evident at the annual scale. Linear scaling with I * and the use of all possible permutations of tracer surveys indicates that flood sequence is not exerting a first-order control on tracer displacement at annual timescales.
The agreement of observed tracer displacement distributions from the Mameyes with models and laboratory data (Lajeunesse et al., 2010) gives us hope that the results here are general. Our earlier observations of tracer dispersion provided an indirect method for inferring particle rest times (Phillips et al., 2013); thus, results from this tracer study may provide the basis for future probabilistic modeling of longterm bed load transport, for which the distributions of particle steps and rests are required input parameters (Zhang et al., 2012). As pointed out above, although larger floods have occurred in the historical record of instantaneous discharge, the associated stresses were no more than 1.28 times the largest observed flood. From the USGS records of annual maximum flood peaks (1967-2012, discontinuous), the largest flood is only 1.35 times our observed largest value of U * . The rapid decay of frequency of occurrence for flows above threshold (Fig. 6a) indicates that transport conditions observed during our study are representative of the river system. Partial bed load transport during near-threshold conditions is also consistent with expectations from equilibrium channel theory for gravel rivers (Parker, 1978;Parker et al., 2007). Indeed, the channel depth inferred from the hydrograph for the "characteristic flood" on the Mameyes -i.e., the peak value in the impulse distribution -agrees with an independent estimate of bankfull flow deduced from vegetation markers, channel morphology, and flow frequency analysis (Pike, 2008;Pike and Scatena, 2010).
We now turn to data from the Bisley 3 tracer deployment, which can serve as a critical test of the generality of the impulse framework and tracer displacement results. Because it is currently unknown how many tracer particles are required to produce accurate statistics, we note here that we only analyze the limited number of Bisley 3 tracers alongside the larger set of Mameyes tracers. Until such a number is known we caution that researchers should not solely rely on a limited set of tracers to inform their results. Given the small number of tracers deployed and the particularly variable stream profile (Fig. 3d), it is intriguing that the Bisley 3 data fall on a well-defined linear relation when plotted Figure 13. Collapse of the tracer data for the first (red +) and second (blue x) Mameyes tracer populations, and Bisley 3 tracers (green o). The black line represents a linear relationship between the mean particle displacement and the dimensionless impulse over the dimensionless friction factor (Eq. 5). The grouping of the date is a result of surveys in three separate years.
against I * (Fig. 10). The offset between the two field sites indicates that, for equivalent values of I * , tracer particles at the Mameyes field site have traveled farther than those in the Bisley 3 stream. This could result from enhanced particle trapping and hiding effects, or greater flow resistance due to the rougher bed. We attempt to collapse the Mameyes and Bisley 3 data onto a single curve by accounting for each of these two effects separately. We assess particle hiding effects using a simple hiding function (Einstein, 1950;Wilcock and Crowe, 2003), which does not produce a collapse of the data. To test the effect of flow resistance, we calculate the dimensionless friction factor k f using a modified Keulegan equation that was found to provide a reasonable fit to a large compilation of field data (Ferguson, 2007): where U is the flow velocity (m s −1 ), and K = 0.41 the von Karman constant. The h r (average from three cross sections) for the Mameyes main channel and Bisley 3 sites is 0.9 and 0.29 m, and the D 84 is 0.31 and 0.55 m, respectively. To reduce variability in the surface grain size counts, the D 84 is determined from an exponential fit to pebble count data collected up and downstream of each study reach. For the Mameyes and Bisley 3 reaches, k f is 0.30 and 9.27, respectively. When values for I * are normalized by √ k f computed using Eq. (5), the mean tracer displacement data for the Mameyes and Bisley 3 streams collapse onto a single curve (Fig. 13). This collapse indicates that, when one accounts for the difference in relative submergence and its effect on flow resistance, the resulting mean particle transport distance is the same. The morphologic characteristics of these two field sites represent end members for the Mameyes watershed, indicating that the linear function X/D = 0.025 I * / √ k f may be a general relationship. This conjecture would be well supported should this relationship be found to hold for tracer studies in other regions. We note here that one can achieve a similar collapse of the data using other recently proposed flow resistance equations as well (see Ferguson, 2007;Rickenmann and Recking, 2011;Ferguson, 2012). A pitfall of the dimensionless impulse is its sensitivity to the determination of U * c . As seen in this manuscript, and in general (Wilcock, 1988), the value of the threshold stress is dependent on the method used to determine it. We have determined two separate values for U * c by using the intercept of the fraction mobile data (U * c1 = 0.22, τ * c1 = 0.023, Fig. 5), and the value that best collapses the long-term mean displacement data (U * c2 = 0.28, τ * c2 = 0.038, Fig. 10). Both of the determined values are within the range reported from field and laboratory data (Buffington and Montgomery, 1997;Mueller et al., 2005;Lamb et al., 2008); however the range of values we have recorded underscores the need for an independent empirical measure of the threshold of motion.

Sediment sorting
The pattern of smaller particles having larger displacements (Fig. 11b) is near universally observed in tracer studies (e.g., Church and Hassan, 1992;Ferguson and Wathen, 1998;Hodge et al., 2011;Scheingross et al., 2013;Schneider et al., 2014). Here we offer a potential explanation for these observations. Results indicate that the modal displacement length for particles of all sizes is approximately one step (Figs. 7a,b,9), and that the dimensionless particle step length depends only weakly on particle size (Fig. 11b). The latter result appears to support laboratory experiments that show that the largest particles travel the farthest for rivers with steep to moderate slopes (Solari and Hill et al., 2010). This effect has been attributed to particle inertia for narrow unimodal sediment size distributions (Solari and , and to bed roughness effects for wider grain size distributions (Hill et al., 2010). The observation that the majority of tracer particles' displacements match the expected step length (Fig. 11b) combined with grain-size-dependent entrainment (Fig. 11a) suggests that the relatively larger displacements for smaller particles result from a greater frequency of entrainment events during a flood. In other words, the largest particles appear to take one step during a flood due to a single entrainment event, while small particles may take multiple steps. Our field results suggest that sorting happens through smaller particles possessing a higher probability of re-entrainment, rather than possessing longer step lengths as compared to larger particles. This may be a consequence of near-threshold transport conditions in the Mameyes, but more work is needed to support this hypothesis.
The cumulative sorting results over annual timescales appear to substantiate aspects of the self-similar sorting theory of Fedele and Paola (2007), and are in general agreement with earlier laboratory experiments (Paola et al., 1992;Paola and Seal, 1995;Seal et al., 1997;Toro-Escobar et al., 2000). The two tracer plumes in the Mameyes had significantly different front positions (Fig. 12a) at the end of the study, as they were emplaced in different years; the first population terminated at X = 1 km, while the second ended at X = 1.2 km. However, the two populations behave dynamically similar when their travel distances are rescaled by their extraction lengths (Fig. 12b-e). The decrease in the standard deviation for tracers is offset from that of the river bed, and results from the narrow grain size distribution of the tracers when compared to the substrate of the river (Fig. 1c). Accordingly, we normalize the tracer and river-bed data by their initial standard deviation at the start of the depositional section of the river (X = 0 km); the result is that tracer sorting appears to track that of the river bed when the initial particle size population and distance traveled are taken into account (Fig. 12b-e). Data show that the initial establishment of the sorting profile can be quite rapid, in that the tracer particles behave dynamically similar to the stream despite different residence times within the river. To our knowledge, this tracer study represents the first active field confirmation of the selective deposition theory (Paola et al., 1992;Paola and Seal, 1995;Seal et al., 1997;Toro-Escobar et al., 2000;Fedele and Paola, 2007).
The full model predicting the mean concentration of the sediment plume downstream consistently underpredicts the data further downstream (Fig. 12e). There are several reasons to expect that a finite population of sediment tracers might not follow the model in Eq. (4). The low values of σ at the leading edge (X * = 1) of the tracer plume may have reached a limit where the size differences between tracers and the stream is negligible and sorting ceases for this mixture. For the stream, the underprediction is potentially from undersampling at the downstream end due to an artificial truncation caused by anthropogenic modification of the stream, which results in a higher value of C 1 in Eq. (3). The scaling exponent is also fairly sensitive to the determination of L in calculating X * in this short system. In the case of the stream bed a small shift in L can steepen or elongate the sorting profile resulting in larger or smaller scaling exponents, respectively. Another factor complicating the scaling of the sorting in the Mameyes River is that the distance from the mountains to the ocean is particularly short (5.95 km), resulting in an abrupt truncation of the sorting profile somewhere within the final kilometer of the river. Given the short time that it took for the tracer particles to sort, and the minimal decline in D 50 downstream, one might expect that there should be a rapidly prograding gravel front (Parker and Cui, 1998); however this front may be arrested due to Holocene sea level rise (Toscano and Macintyre, 2003).

Conclusions
In this paper we have presented field results on bed load tracer displacement data at the event to annual timescales, and used simple theory to rationalize the displacement scaling of tracer particles, show how a heterogeneous population of particles sorts downstream, and explore the implications of these findings on the statistical scaling of the hydrograph. At the scale of single floods, the distribution of particle displacements is well described by an exponential distribution. Close agreement with laboratory data and theory (Lajeunesse et al., 2010) suggests that these displacements represent the scaling of the fundamental particle step length. We infer that, for near-threshold floods, the most probable transport distance is one step length. Cumulative displacement over many floods reinforces this finding, with data showing that tracers remain in the partial transport regime for a range of flow conditions. We test the applicability of the impulse framework using data from two streams of very different morphologies, and find that tracer displacement data collapse onto a single linear relationship when differences in critical Shields stress and flow resistance are accounted for. For particle sorting, we find that downstream fining emerges after a series of floods. Sorting seems to result from a slight difference in size-dependent particle entrainment at the flood scale. We find that tracers sort to the limit of sorting present in the stream bed. Both the tracers and the stream bed have the same scaling when accounting for the distance each has traveled, as well as the initial statistics of the tracer and stream grain size distributions, respectively. These observations serve as an active field validation of the selective deposition sorting model (Paola et al., 1992;Fedele and Paola, 2007). Finally, we show that the magnitude-frequency distribution of flood stress in the Mameyes is exponential for flows exceeding the threshold of motion. The average stress for flows exceeding critical is approximately 1.4 times the critical Shields stress, and represents the stress of maximum geomorphic work. In addition, the distribution of dimensionless impulse has a welldefined peak coincident with the flood of maximum geomorphic work, indicating that the channel is adjusted to a characteristic flood impulse. We believe that tracer dynamics observed in the Mameyes River are characteristic of many gravel rivers, because many gravel streams are adjusted such that bankfull floods exert a stress that is only slightly in excess of the threshold for entrainment (Andrews, 1984;Pitlick and Cress, 2002;Torizzo and Pitlick, 2004;Mueller et al., 2005;Parker et al., 2007). A caveat for all of our results, however, is that caution should be exercised when considering reported numerical values due to the difficulty in independently determining the threshold of particle motion. We emphasize that this remains one of the most critical problems in determining coarse-grained sediment mechanics in natural rivers.