Percolation observations in an arid piedmont watershed and linkages to historical conditions in the Chihuahuan Desert

. A critical hydrologic process in arid and semiarid regions is the interaction between ephemeral channels and groundwater aquifers. Generally, it has been found that ephemeral channels contribute to groundwater recharge when stream ﬂ ow in ﬁ ltrates into the sandy bottoms of channels. This process has traditionally been studied in channels that drain large areas (tens to hundreds of square kilometers). Since the water table in arid and semiarid regions is typically far from the surface, measured stream ﬂ ow losses or percolation into the deep vadose zone is equated to groundwater recharge. In this study, we use a water balance approach to estimate deep percolation in a ﬁ rst-order, instrumented watershed (4.7 ha) on a pied-mont slope of the Jornada Experimental Range (JER) in the Chihuahuan Desert. Results indicate that runoff generated within the piedmont slope contributes signi ﬁ cantly to deep percolation. During the short-term 6-yr study period, we estimated 385 mm of total percolation, 62 mm/yr, or a ratio of percolation to rainfall of 0.26. Based on the instrument network, we identi ﬁ ed that percolation occurs inside channel areas when these receive overland sheet ﬂ ow from hillslopes. We observed less stream ﬂ ow leaving the watershed as compared to percolation during the study period, leading to an outlet stream ﬂ ow to rainfall ratio of 0.02. Using long-term data sets available from the JER, we estimate that over the last 100 yr, 48 mm/yr of percolation occurs at the study site, a ratio of percolation to rainfall of 0.19. We then scale this up to determine the contribution of similarly sized watersheds on the piedmont slopes of the JER. These observations high-light the importance of arid piedmont slopes for generating groundwater recharge, in particular during above-average rainfall periods.


INTRODUCTION
Recharge in arid and semiarid regions represents a challenging topic for the hydrologic sciences (Phillips et al. 2004).Due to its practical importance for the sustainable management of groundwater aquifers, there has been a concerted effort in recent years to understand this process.The traditional view of recharge in arid regions has been as a sum of several different pathways, including mountain block recharge (Mailloux et al. 1999), mountain front recharge (Manning and Solomon 2003), diffuse recharge (Stephens 1994, Small 2005), and ephemeral channel recharge via channel transmission losses (Phillips et al. 2004).More recent studies combine mountain block and mountain front recharge into a single term (Ajami et al. 2012, Meixner et al. 2015).Additionally, the important role of vegetation in moderating, and often negating, recharge has been recognized.Field studies employing lysimetry (Gee et al. 1994, Levitt et al. 1999), time domain reflectometry (Smith et al. 1995), and onsite thermocouple psychrometry (Andraski 1997, Scanlon et al. 1999) in naturally vegetated environments report no deep percolation (or diffuse recharge) beneath the top 1-2 m of soil due to vegetation uptake.Gee et al. (1994) found that when vegetation was removed, diffuse recharge rates increased by >50% of the annual rainfall on a piedmont slope in the Jornada Experimental Range (JER).
Recharge in the southwestern United States is dominated geomorphically by the basin and range physiographic province where alluvial fans, piedmont (bajada) slopes, and closed basins are common (Monger et al. 2006).It is recognized that topographic lows, such as playa lakes and ephemeral channels, likely contribute significantly to recharge (Phillips et al. 2004).Many studies have focused on the role played by playas in facilitating recharge due to their periodic inundation (Wood and Sanford 1995, Scanlon and Goldsmith 1997, Yechieli and Wood 2002, McKenna 2016).Traditionally, the role played by channels in controlling recharge has been quantified by measuring channel transmission losses between two streamflow gauging stations (Renard et al. 1993, Goodrich et al. 1997, 2004).Other approaches also exist.For example, Bailey (2002) estimated ephemeral channel losses using temperature monitoring, finding that vertical transport during streamflow is most sensitive to the channel hydraulic conductivity.This was expanded upon by Blasch et al. (2006) who showed that infiltration rates are high (~1000 m/ d) at the onset of streamflow in ephemeral channels, but decay by several orders of magnitude within the first hour of channel transmission losses.Additionally, soil water content measurements are a common method for quantifying recharge in ephemeral channels (Dowman et al. 2003, Shanafield andCook 2014), but examples of this approach tend to focus on large channels that transport significant amounts of runoff generated in upstream areas.To our knowledge, similar efforts to quantify groundwater recharge in first-order watersheds that primarily generate runoff locally have yet to be conducted in arid and semiarid regions.
A tool that has been underutilized in estimating ephemeral channel losses is the water balance approach, mainly due to the necessity of highresolution measurements of simultaneous processes in a watershed.In addition to the data requirements, another limitation of the water balance approach for calculating deep percolation (P) is that as a residual estimation, it requires differencing two large fluxes-rainfall (R) and evapotranspiration (ET)-subject to measurement and spatial aggregation errors.Many studies in arid and semiarid regions have found that the local (soil) water balance is dominated by ET such that ET/R % 1, in particular for flat surfaces with limited runoff generation potential (Dugas et al. 1996, Scott et al. 2004, 2006, Scott 2010, Glenn et al. 2015).However, prior efforts in shrubland sites with a large amount of bare soil (>60%) have quantified lower ET than in other ecosystems (Dugas et al. 1996) and reported values of ET/R < 1 allowing for the generation of streamflow and/or deep percolation (Glenn et al. 2015).When streamflow or deep percolation is sufficiently large, then the sensitivity of the water balance approach to errors in the measurement of R and ET is expected to be smaller.
For instance, water balance studies have been conducted in the semiarid Walnut Gulch basin (140 km 2 ) in southeastern Arizona, which contains Chihuahuan Desert shrublands (Kustas and Goodrich 1994, Scott et al. 2004, 2006).Using long-term observations, Renard et al. (2008) quantified the annual water balance in the watershed, indicating that ET accounts for 93% of R, channel transmission losses represent 6% of R, and streamflow at the outlet (Q) accounts for only 1% of R. In this conceptual model, transmission losses of onsite runoff occur mainly in channels with large upstream areas (Renard et al. 2008), determined by differencing streamflow volumes from paired gauging stations.In addition, it has been recognized that ephemeral channel recharge can be a more significant component of the water balance during wet years (Goodrich et al. 2004).Nevertheless, the contribution of first-order watersheds on the arid piedmont slope to groundwater recharge has not been determined, in part due to a lack of water balance studies in smaller areas.Another factor is that the geological conditions of small watersheds and their internal channels can vary significantly in terms of their connectedness to groundwater aquifers.As a result, channel transmission losses (deep percolation) often need to be taken as a proxy for groundwater recharge (Wilson and Guan 2004) despite a lack of direct measurement of the flux into the deep water table .This assumption is likely more valid for ephemeral channels lacking a phreatophytic vegetation community along banks that would consume transmission losses.
In this study, we aim to quantify deep percolation in a small watershed located on an arid piedmont slope in the Chihuahuan Desert using the water balance approach applied to measurements from a dense instrument network over the period of 1 July 2010 to 30 September 2016.We obtain deep percolation (P) as a water balance residual.For the purposes of this study, we define P as the exchange of water between the watershed soil surface (40 cm depth) and deeper subsurface materials.Positive (negative) values of percolation indicate a downward (upward) water flux.The soil depth was selected due to the common presence of a caliche (or calcium carbonate) layer serving as a semi-impermeable lower boundary and limiting the placement of soil moisture sensors (Templeton et al. 2014).Through the application of the water balance, we derive monthly and annual estimates of deep percolation during the study period and assess their sensitivity to measurement errors.Using the instrument network, we then explore the potential mechanisms leading to deep percolation and describe their occurrence within the watershed.By linking the deep percolation to long-term rainfall data, we then extend the temporal context of the study to assess how wet years might influence groundwater recharge.We also identify a number of sites on piedmont slopes of the JER where deep percolation might play a significant role.Given the location of these first-order watersheds, deep percolation estimates should be considered as piedmont slope recharge, a flux that has not been commonly quantified (Phillips et al. 2004) or included in conceptual models of arid watersheds (Renard et al. 2008).

Study site
The JER, located ~20 km north of Las Cruces, New Mexico, USA, within the Chihuahuan Desert, has been used continuously as an experimental rangeland since 1912.The basin and range province in the region consists of north-south trending mountain ranges and broad intervening desert basins that are the product of the Rio Grande rift tectonic system (Monger et al. 2006).The major physiographic components of the JER are the San Andres mountains, piedmont slopes, basin floors, and the Rio Grande valley, with elevations that range from 1180 to 2747 m.The piedmont slopes within the JER are composed of alluvial deposits that can be >1200 m thick and have thick unsaturated zones where the depth to groundwater ranges from 90 to 105 m (King and Hawley 1975).This geomorphic setting leads to runoff produced in mountains that is transported via ephemeral channels across piedmont slopes and onto the basin floor where playas may form.Large drainage systems spanning from mountains to playas also interact with smaller channels that form directly in the piedmont slopes.An example of these small ephemeral channels is shown in Fig. 1, characterized by a length of 400 m, width of 0.5 m near the outlet, and mean slope of 3.5%.The area draining to the ephemeral channel is a first-order (4.7 ha) watershed on the piedmont slope containing soils derived from the San Andres mountains that are primarily sandy loam in texture with high gravel contents (Anderson and Vivoni 2016).Distinct hillslopes facing north, south, and west drain an upper flat section toward the westward flowing channel.On the hillslopes, a well-defined caliche layer can be found at a depth of 30-40 cm (Templeton et al. 2014), while the channel sediments are deeper than hillslope soils (typically >100 cm).
The mean annual rainfall in the watershed (2010-2016) of 257 mm/yr is dominated by summer storms (70% of R occurring from July to September).While this average is consistent with long-term rainfall monitoring at the JER Headquarters (1915-2016) and Taylor Well (1980-2016) of 245 and 246 mm/yr, the percentage of annual rainfall in the summer is closer to 53% for the longterm data, as shown in Table 1 (JER Data Catalog: jornada.nmsu.edu/data-catalogs/long-term).Given the rainfall amounts and average annual temperature of 15.3°C, climate conditions are cool and arid (K€ oppen zone BWk).Based on the Thornthwaite classification, the region is arid (zone E) and only rarely crosses the threshold into semiarid (Wainwright 2006), consistent with a low ratio of rainfall to potential evapotranspiration (R/PET = 0.11).Climate conditions and the process of woody plant encroachment in the last 150 yr (Gibbens et al. 2005) have led to a mixed shrubland ecosystem at the site consisting of creosote bush (Larrea tridentata), honey mesquite (Prosopis glandulosa Torr.), other shrubs (Parthenium incanum, Flourensia cernua, and Gutierrezia sarothrae), and several grass species (Muhlenbergia porteri, Pleuraphis mutica, and Sporobolus cryptandrus).Shrubs at the study site are typically small (~1 m in height) and widely spaced, leading to a high bare soil coverage of 65% (Fig. 1).To quantify these conditions, a high-resolution (1 m) image analysis of terrain and vegetation properties by Templeton et al. (2014) was used to derive hillslope and channel network properties and plant species cover in the watershed (Table 1).2016) describe monitoring efforts that began in 2010 at the watershed with the establishment of a dense network of rainfall, runoff, soil moisture, and soil temperature sensors, as well as meteorological, radiation, and energy flux measurements at an eddy covariance (EC) tower (Fig. 1).This brief description is focused on the measurements and data processing methods used to calculate the watershed water balance.Rainfall (R) was measured using up to four tipping-bucket rain gauges (TE525MM; Texas Electronics, Dallas, Texas, USA) to construct a 30-min resolution spatial average based on Thiessen polygons.Streamflow (Q) at 1-min and 30-min intervals was measured at the watershed outlet using a Santa Rita supercritical runoff flume (Outlet Flume; Smith et al. 1981), a pressure transducer (CS450; Campbell Scientific, Logan, Utah, USA), and an in situ calibration (Turnbull et al. 2013).Additionally, channel runoff was obtained at three interior locations using smaller flumes (Wainwright et al. 2002), pressure transducers (CS450; Campbell Scientific), and an in situ calibration (Templeton 2011).Evapotranspiration was obtained at 30-min resolution using the EC technique that employs a three-dimensional sonic anemometer (CSAT3; Campbell Scientific) and an open-path infrared gas analyzer (LI-7500; LI-COR, Lincoln, Nebraska, USA) installed at a 7 m height.Surface energy fluxes were calculated following the procedures described in Templeton et al. (2014) and Anderson and Vivoni (2016).Watershed average soil moisture and temperature measurements at 30-min resolution were obtained using soil dielectric probes (Hydra Probe; Stevens Water, Portland, Oregon, USA) organized as profil e s( s e n s o r sp l a c e da t5 , 15, and 30 cm depths) in three transects along each major hillslope (Fig. 1b).At the EC tower, soil moisture sensors were installed at 5, 15, 30, and 50 cm depths, as the deeper soils on the upper flat section allowed installation of a deeper sensor.Additionally, a cosmic-ray neutron sensor (CRNS; H y d r o i n n o v a ,A l b e q u e r q u e ,N e wM e x i c o ,U S A ) was installed in 2013 to obtain a spatially averaged soil moisture value (h CRNS ).Data from the soil moisture probes were averaged by depth and location to obtain a spatially averaged quantity (h SN ) consistent with the CRNS measurement (Schreiner-McGraw et al. 2016).CRNS data were used to gap-fill the soil moisture observations when necessary.To complement these long-term measurements, a profile of soil moisture sensors was deployed in the main channel on 4 August 2016, at depths of 15, 30, 50, and 100 cm (Fig. 1c).Notes: Vegetation species have been generalized into three types.Rainfall values are long-term averages  with seasons defined as fall (October-December), winter (January-March), spring (April-June), and summer (July-September).

Calculating percolation through independent methods
The watershed water balance can be expressed as: where Dh is the change in volumetric soil moisture over the time interval Dt, Z r is the soil depth, R is rainfall, ET is evapotranspiration, Q is streamflow at the outlet, and P is deep percolation below Z r .All terms in Eq. 1 are expressed in mm per Dt as spatially averaged values over the watershed.By measuring Dh,R,ET, and Q at high temporal and spatial resolution, deep percolation can be obtained as a residual from (Eq. 1).Positive P values indicate water is percolating from the shallow soil surface to the deeper vadose zone, while negative values suggest that plants are transpiring water from beneath Z r ,d e fined here as 40 cm, as described previously.We evaluated the watershed water balance at a monthly time step (Dt = 1 month) to ensure that rainfall pulses are redistributed within the watershed or evapotranspired within the shallow soil depth.A total of 4 months in the 75-month record (2010-2016) experienced instrument failure and were excluded from the analyses.The water balance was also evaluated at an annual scale (and for incomplete years), and cumulative percolation was obtained over the entire record.
In order to calculate the monthly or annual water balance components, high-temporal (30min) resolution observations from the instrument network were aggregated first to a daily scale.Fig. 2 presents the spatially averaged daily values of the terms of (Eq. 1) used to estimate deep percolation over the study period, including a comparison of soil moisture measurements from the sensor probe network and the CRNS sensor during the overlapping period.Notice how the terms of (Eq. 1) are dominated by summer storm events from July to September that are part of the North American monsoon system (Douglas et al. 1993, Vivoni et al. 2008), with the majority of rainfall (71%), ET (53%), and streamflow (100%) in this period.Summer periods in 2013 and 2014 had particularly large rainfall amounts, exhibiting high soil moisture values that matched well in the h SN and h CRNS methods (Schreiner-McGraw et al. 2016).Fig. 2 is also illustrative of the difficulties of evaluating the watershed water balance at a daily temporal resolution to estimate deep percolation, namely: (1) R is infrequent, with some large events having >50 mm/d, but many days characterized by no rainfall, (2) ET is a nearly continuous process occurring at low rates (0 to ~3 mm/d), but varies strongly with the amount of soil moisture, (3) Q at the outlet is a rare occurrence in response to specific rainfall and soil moisture conditions, and (4) h exhibits variability stemming from the characteristics of individual storm events as well as the seasonality inherent in the rainfall and radiative forcing to the land surface.
We accounted for potential errors in the measurement of the major terms of (Eq.1), R,E T , and Q, to quantify the uncertainties associated with the percolation estimate, following Scott (2010).For R, the effect of gauge undercatch has been estimated as 5% of the measured rainfall (or 1.05R; Duchon and Essenberg 2001) due to the effects of high winds and large rainfall intensities on rain gauge collection (Sieck et al. 2007).For ET, we applied the Bowen ratio (B = H/kET) method of Twine et al. (2000) to force energy balance closure at the daily scale.The closure error was found to be 19% during the study period (Schreiner-McGraw et al. 2016) and is consistent with reports in other regions (Wilson et al. 2002, Scott 2010).By assigning the error into latent (kET) or sensible heat (H) flux components based on a seasonally averaged B, a closure-forced ET (ET f ) of a higher magnitude was derived.For Q, we considered a worst-case scenario of 100% error (or 2Q) to provide a conservative estimate of deep percolation.Different combinations of the water balance measurement errors were applied across all years to obtain a range of possible annual values of P and its cumulative amount over the entire period.
We compared percolation estimates derived from the watershed water balance to soil moisture measurements in the ephemeral channel obtained over a short period with available data (4 August to 30 September 2016).Streambed percolation (P SB ; mm) was obtained following Shanafield and Cook (2014) using soil moisture measurements at 15, 30, 50, and 100 cm depths: where Z c is the depth of channel sediments (assumed to be 100 cm) and h c is the depthaveraged channel soil moisture, while W c and L c are the average width and total length of the channel, respectively, and A t is the total area of the watershed.Additionally, the channel soil moisture measurements were compared to corresponding observations at the EC tower (15, 30, and 50 cm depths) to assess the movement of infiltration fronts into channel and upper flat sites in the watershed.Similar analyses were conducted with sensors from the three transects to characterize the occurrence of vertical and lateral water movement within hillslope soils.Finally, channel transmission losses were estimated using observations at the three internal runoff flumes and the outlet flume via differencing of storm-based runoff volumes (Goodrich et al. 1997).

Assessing historical and spatial implications
We attempt to place the study site observations in a broader spatiotemporal context to assess the importance of deep percolation in arid piedmont slopes.For this purpose, we construct a relationship between observed monthly P and R, as performed in prior studies (Ajami et al. 2012), and fit a piecewise linear regression using a nonlinear optimization approach.By assuming the relation holds over longer time periods, we generate deep percolation estimates from monthly rainfall data available at JER Headquarters  and Taylor Well (1980Well ( -2016) ) to identify periods of above-and below-average P in the historical record.Additionally, we utilize a 10-m digital elevation model of the JER to delineate all first-order watersheds on piedmont slopes with a range of areas matching the study site (4-5 ha).Under the assumption of similar percolation values, we estimate the potential contribution of first-order watersheds to piedmont slope recharge and compare these estimates to other recharge sources.

RESULTS AND DISCUSSION
Percolation estimates and their uncertainty Fig. 3 shows monthly estimates of deep percolation obtained as a water balance residual (accounting for rainfall undercatch as 1.05R and closure-forced ET as ET f ) as well as the cumulative P during the study period.Annual P estimates are also indicated in relation to the rainfall during each period.For years with R less than average (R < 257 mm/yr), monthly P generally fluctuates between À20 and +20 mm, depending on amounts of monthly R and ET.Small values of negative monthly P indicate plant consumption (ET) of water stored below 40 cm and typically occur during the winter and spring (January-June).Positive monthly P occurs during the summer monsoon (July-September) when large rainfall events exceed the capacity of plants to consume water via ET and the streamflow (Q) from the watershed within a month is small.For the years 2013 and 2014 with aboveaverage rainfall, monthly P exceeded +20 mm during five months, reaching values of 78 and 95 mm during September of each year.A concentration of high rainfall amounts during the monsoon promoted percolation.For instance, three storms in summer 2014 (43, 62, and 51 mm of total rainfall over an average duration of 11 h) were responsible for 52% of the annual R.A sa result, these two years accounted for 75% of the cumulative P of 385 mm and had a ratio of P/R of 0.45, significantly larger than the overall ratio during the study period of P/R = 0.26.In addition, P/R in 2013 and 2014 was substantially larger than in 2015 (P/R = 0.18) with nearly average rainfall (R = 235 mm) but a more even distribution during the year (34% of R from July to September; Fig. 2).
Table 2 presents annual measurements of each water balance component and compares percolation estimates under different measurement error assumptions.Evapotranspiration is the largest component of R with a total ET of 1050 mm over the study period (ET/R = 0.72), while the use of a closure error correction increases this  value to ET f = 1158 mm or ET f /R = 0.80.ET and ET f show a slight increase in 2015 and 2016, despite the lower-than-average rainfall in those periods, suggesting a small amount of consumption of stored water.In contrast, streamflow (Q = 39 mm, Q/R = 0.02) and soil moisture changes (Z r Dh/Dt = À54 mm, Z r Dh/RDt = À0.04)have a small contribution to the water balance during the study period.Moreover, large differences in R and ET during above-average rainfall years (2013 and 2014) are not reflected in substantial increases in Q or Z r Dh/Dt.These comparisons suggest that for summer monsoons when rainfall is in a large excess of concurrent ET losses, there is a significant amount of deep percolation.The different combinations of measurement error assumptions shown in Table 2 do not appreciably alter this outcome.For instance, when all losses in the water balance are maximized (P=RÀ ET f À 2Q À Z r Dh/Dt, no undercatch, closure-forced ET, and 100% error in Q), deep percolation is still 19% of total rainfall (or 7% lower than the standard estimates; bottom row in Table 2).Overall, measurement errors lead to a range of total P from 273 to 493 mm or P/R of 19-34%, with the majority of the percolation occurring during summer periods with the above-average rainfall.

Evaluation of percolation mechanisms
Given the large estimates of deep percolation, Fig. 4 uses the dense instrument network to explore relations between P and other hydrologic quantities that might explain the underlying mechanism.A strong positive relation of a power-law form (P = 1.32Q 2.23 , p = 0.02, RMSE = 27.5 mm) was found between monthly P and streamflow at the outlet for months when Q > 0 (Fig. 4a).Thus, for summer months when streamflow occurs (n R = 15 months), positive residuals of the watershed water balance (P > 0; Fig. 3) are correlated with streamflow yield.This suggests that deep percolation and streamflow are both linked to hydrologic processes occurring in the channel network.In contrast, no significant relation was obtained between P and monthly soil moisture changes (Z r Dh/Dt)f o rm o n t h sw h e nQ > 0, suggesting that positive residuals (P > 0) are not related to water storage changes in the hillslopes.Similarly, we found no significant relation between P and Z r Dh/Dt for months without streamflow (n NR = 59 months), where both positive and negative water balance residuals are possible (Fig. 4b).This is consistent with the low sensitivity of annual values of Z r Dh/Dt to the year-to-year variability in rainfall (Table 2).As a result, the water balance residuals do not appear to be related to water storage changes in hillslope soils during months with large P (average over n R of 30 mm/month) when rainfall is in excess of ET and leads to streamflow, or across other periods with small positive and negative P (average over n NR of 0 mm/month).As noted previously (Fig. 2), negative P values indicate an upward movement of water from depths greater than Z r which could be associated with ET at the monthly scale.For those months with Q = 0, we found a significant negative relation between P andET(P = À0.24ETÀ 0.88; p = 0.02).For those months with P < 0andQ = 0, the negative relation between P and ET is more statistically significant (P = À0.29ETÀ 1.77; p = 0.0001).This suggests that negative values of P are related to the consumption of deeper sources of water by plants that negate diffuse recharge on hillslope soils, while large positive values of P are related to direct recharge of the deeper subsurface within the permeable channel sediments.
To illustrate differences between channel and hillslope conditions, Fig. 5 presents soil moisture (h) observations at three depths (15, 30, and 50 cm) at the EC tower (hillslope) and the main channel (including an additional 100 cm sensor), along with rainfall and streamflow data.A storm event on 5 August 2016 (16 mm over 5 h) resulted in streamflow at the outlet of 3.7 mm with a peak 23 min after rainfall started.While small storms occurred the previous week (total rainfall of 4 mm), overall conditions in July were drier than average (13% of average R) as noted in the soil moisture data.On hillslope soils, limited infiltration was observed at 15 cm depth in response to the storm, and no change in soil moisture was noted at 30 or 50 cm depths.In sharp contrast, the channel sensors responded to the runoff generated during the storm event at all depths as well as to a number of subsequent rainfall periods that did not generate outlet streamflow.Furthermore, deeper sensors showed a slower soil moisture recession than shallower ones, indicating vertical water movement through the channel sediments over time.Percolation in the channel, even in the absence of outlet streamflow, suggests that topographic lows in first-order watersheds are locations where runoff losses occur and are responsible for positive water balance residuals (P > 0) when rainfall inputs exceed ET.This was confirmed by the good match between streambed percolation (P SB = 10 and 0 mm) and the deep percolation estimate (P = 11.9 and 0 mm) in the period of August-September 2016.
To confirm that hillslope soils exhibit limited infiltration, soil moisture observations from the sensor network were inspected.Analysis of soil moisture variations along each hillslope transect following Guti errez-Jurado et al. ( 2007) did not reveal subsurface downslope movement following large storms.Fig. 6 provides an example for a 44-mm storm on 13 September 2013 responsible for a large positive water balance residual (P > 0) inferred to be deep percolation.Over the course of a 5.5-d period, the spatially averaged soil moisture in the sensor network (h SN ) increased in response to the storm from 0.11 to 0.24 m 3 /m 3 at the peak value.Most of the increases in soil moisture, however, were limited to the top 30 cm of soil.Soil moisture sensors at the EC tower showed no increases below this depth.Furthermore, the depth-averaged soil moisture exhibited a nearly uniform spatial distribution during each of the selected time periods.Clearly, the uniformly dry watershed (antecedent condition) became uniformly wet up to a shallow depth (peak soil moisture), followed by a progressive drying that occurs nearly at the same rate within a single transect.Spatial differences, however, are noted among transects for the last period (120 h after peak) when h SN = 0.13 m 3 /m 3 , due to variations in microtopography.Note that the westfacing transect is located in a slight topographic low that promotes higher soil water accumulation.The uniform wetting and drying processes in hillslope transects indicate that subsurface water movement is not a significant mechanism for producing runoff.Furthermore, the generally shallow wetting fronts and drier soil layers at depth suggest that hillslope runoff generation is primarily from an infiltration-excess mechanism that produces overland transport (sheetflow) to the receiving channel reaches at the bottom of the hillslopes.
To verify that overland flow is lost in channel reaches, we inspected measurements at the outlet flume and three internal flumes (Flumes 1, 2, and 3; Fig. 1).Fig. 7 describes the spatial and temporal variations of runoff to infer transmissions losses during individual storms and the entire study period.High-resolution measurements (Fig. 7a) for a storm on 4 September 2014 illustrate that runoff begins and ends earlier in tributaries (Flumes 1 and 3) as compared to the main reach (Flume 2 and Outlet) consistent with having a source from upstream hillslopes.The peak recorded at Flume 2 takes 26 min to arrive at the outlet and runoff lasts for 36 min after the outlet streamflow ceases, indicating transmission losses in the main reach.This is confirmed by inspecting the relation between upstream runoff (sum of Flumes 2 and 3) and outlet streamflow normalized by drainage area (Fig. 7b).When inspected over the events with data available at all flumes, upstream runoff is larger than outlet streamflow by 108%.When the reach between Flumes 1 and 2 is also considered, transmission losses are 21.7 mm over the study period or an average of 4 mm/yr, a magnitude comparable to the measured Q (Table 2).However, this estimate is considered a lower bound on the P estimated as a water balance residual since there are large soil moisture increases within channel reaches prior to the generation and transport of runoff (Fig. 5).Furthermore, a piecewise linear relation is noted between total runoff volume and drainage area (Fig. 7c), with increases in runoff up to about 2.5 ha, after which a slight reduction occurs.This suggests that a threshold exists in drainage area after which transmission losses become important, for instance between Flume 2 and the Outlet, due to the channel size and its storage capacity.The relation also shows that a minimum drainage area of 0.08 ha exists for runoff initiation.This was verified by comparing the mapped channel network based on GPS readings and the channel reaches derived from a 1-m digital elevation model assuming a channel area threshold of 0.08 ha (Fig. 7c inset).

Conceptual model of percolation mechanism
Observational evidence from the dense instrument network suggests that positive water balance residuals are due to deep percolation in channel reaches.Fig. 8 presents a conceptual model for the percolation mechanism on arid piedmont slopes under the influence of small and large rainfall events during the monsoon (Gebremichael et al. 2007, Villarreal et al. 2016).Rainfall events on hillslope soils with limited infiltration capacity (Rossi et al., unpublished manuscript) result in the production of overland transport (sheetflow) via the infiltration-excess mechanism.During large storms when infiltration fronts penetrate to deeper layers (Schreiner-McGraw et al. 2016), a semi-impermeable calcium carbonate layer limits the hillslope storage capacity (Duniway et al. 2010).Sheetflow production on hillslopes occurs in response to high amounts of bare soil (Tromble 1988), the dense gravel cover on hillslope soils (Templeton et al. 2014), and the hydrophobic compounds derived from desert shrubs that coat soil surfaces (Adams et al. 1970, Estell et al. 1994).When sheetflow accumulates in channels, the coarser-textured sediments allow for rapid infiltration of overland flow which fills pore spaces within the channel storage (S c ).For small events, all hillslope runoff can be effectively stored in channel sediments, leading to positive water balance residuals (P > 0) and no outlet streamflow (Q = 0).Small events are responsible for increases in soil moisture in channels and redistribution into the deeper vadose zone, though limited runoff production.For large rainfall events, sheetflow inputs into channels can overwhelm the channel storage capacity, leading to the transport of internal runoff to the outlet as streamflow and the occurrence of transmission losses as transport occurs.For this case, both positive water balance residuals (P > 0) and outlet streamflow (Q > 0) are present.As a result, channel losses in the first-order watershed during small and large events are responsible for the large estimated amounts of deep percolation which are subsequently redistributed in the subsurface and potentially available as groundwater recharge.By identifying piedmont slope recharge as an important mechanism, the conceptual model of Renard et al. (2008) for the water balance of arid and semiarid watersheds can be refined to include this process and explain their term labeled as onsite runoff in more detail.In the Renard et al. (2008) model, onsite runoff is sheetflow that becomes streamflow in first-order watersheds.Results from this study indicate that an additional term onsite percolation should be added to represent sheetflow that percolates before becoming streamflow as illustrated in Fig. 8a.The observational evidence presented here suggests that arid piedmont slopes produce substantial amounts of sheetflow that are stored in small ephemeral channels and serve as direct recharge rather than contributing directly to ET.

Historical and spatial estimates of percolation
We utilize the hydrologic observations in the watershed to extend the temporal and spatial context of the deep percolation mechanism in arid piedmont slopes.Fig. 9 presents a relation between monthly rainfall and percolation derived from the dense instrument network.A nonlinear optimization approach (Gill et al. 1981) was used to construct a piecewise linear regression with an identified breakpoint at R = 62.8 mm/month.For months with higher  rainfall, deep percolation increases linearly in a positive fashion (P = R À 44.33) with a slope of 1, suggesting that once the storage capacity of the watershed has been exceeded, all excess rainfall becomes percolation.For months with rainfall lower than the threshold, deep percolation can assume low values near zero (P = 0.37R À 4.76).This is consistent with prior analyses that showed that P fluctuates between À20 and +20 mm for most months, while summer periods with high rainfall allow for large increases in P. A sensitivity analysis to the selection of the threshold (from 40 to 70 mm/month) revealed that the ratio of P/R = 0.26 in the watershed did not change appreciably during the study period.This approach is consistent with linear relations often used to model mountain recharge in arid and semiarid regions as a function of annual rainfall (Anderson et al. 1992).For instance, Ajami et al. (2012) used a recharge relation with a summer season threshold of 130 mm for mountain areas in southern Arizona.Moreover, the watershed instrument network allowed for constructing a rainfall-recharge relation at a high (monthly) resolution that is specific to arid piedmont slopes.
The monthly relation between rainfall and deep percolation derived from the first-order watershed is applied to long-term data from two additional rain gauges at the JER to derive annual estimates of piedmont slope recharge.Fig. 10a shows the temporal variability of annual P based on the three rain gauge records, resulting in P/R of 0.19 or P = 48 mm/yr for the long-term sites (100 yr at JER Headquarters and 30 yr at Taylor Well).Notably, the presence of wet years can lead to significant amounts of deep percolation, with above-average P occurring 63% of the time.Wet years at the watershed site, such as 2013 and 2014, are also common in the long-term records, indicative of the frequent presence of high rainfall during certain summer months that promote large amounts of channel losses.In addition, a significant positive trend of 0.43 mm/yr (p = 0.02) is present in annual deep percolation from 1915 to 2015.Furthermore, annual P < 0 is infrequent and three or more consecutive years with negative P do not occur in the long-term record.This indicates that annual water losses to the deeper subsurface through channel sediments are not depleted by deep-rooted vegetation in subsequent years, as suggested to occur for diffuse recharge on flat terrain (Gee et al. 1994, Scanlon et al. 1999).As a result, piedmont slope recharge driven by above-average summer rainfall is likely to largely remain within the subsurface during subsequent dry periods.We have shown that the firstorder watershed likely contributes significantly to groundwater recharge.If this is a common feature, the contribution of first-order watersheds to aquifer recharge would be potentially widespread since the site is not unique, as depicted in Fig. 10b which shows similarly sized (4-5 ha) watersheds located on piedmont slopes in the JER (a total number of 1,737 watersheds or 2.2 sites per km 2 ).Assuming an annual averaged P of 48 mm/yr at all locations yields an annual potential recharge to the groundwater system of 3.8 Mm 3 /yr, which is significantly larger than recent estimates of playa recharge on the basin floor (65,000 m 3 /yr) based on the playa area and annual recharge values of 5 mm/yr (McKenna 2016).This suggests that firstorder watersheds on piedmont slopes are important contributors to groundwater recharge in arid and semiarid regions.

CONCLUSIONS
Based on multiple lines of evidence, we identified that first-order watersheds have the potential to contribute significantly to groundwater recharge on arid piedmont slopes at the JER.The water balance approach using a dense instrument network afforded the ability to identify significant losses during periods when rainfall was in large excess of ET and streamflow yield.Spatiotemporal analyses of hillslope soil moisture, channel percolation, and channel runoff were useful in testing different percolation mechanisms.Future studies should directly measure percolation at depths >1 m at various locations in the channel network.A new conceptual model for piedmont slope recharge emerged whereby the limited infiltration capacity of hillslope soils promotes overland sheetflow that is stored in small receiving channel reaches leading to deep subsurface storage.Under high antecedent wetness or a large rainfall event, the channel storage capacity is exceeded and both deep percolation and streamflow yield are noted.The concept of piedmont slope recharge connects previously identified mechanisms of diffuse recharge occurring in upland soils and direct recharge within large downstream channel systems, thus enriching knowledge in arid and semiarid regions.We analyzed historical records to show that percolation is a common occurrence due to the number of wet years.While the study period captured wet years well, it may not have been varied enough to depict the ET response during prolonged dry periods.Nevertheless, the common nature of firstorder watersheds on piedmont slopes and the long-term variability in rainfall promoting deep percolation suggest that this groundwater recharge mechanism is likely significant for the overall water budget of landscapes in the basin and range province of the southwestern United States.

Fig. 1 .
Fig. 1.(a) Study site location.(b) One-meter aerial photograph with the watershed boundary and stream network as well as the sensor locations.Photographs of (c) the main channel, (d) the south-facing hillslope where sheetflow commonly occurs, and (e) the upper flat surface near eddy covariance tower.

Fig. 2 .
Fig. 2. (a) Daily rainfall, streamflow, and evapotranspiration.(b) Daily spatially averaged soil moisture from the sensor network (h SN ) and the cosmic-ray neutron sensor (CRNS) method (h CRNS ).Gaps in the data indicate periods of instrument malfunction or data collection problems, and the shaded regions indicate the North American monsoon (July-September).

Fig. 3 .
Fig. 3. Monthly percolation (bars) and cumulative percolation (line) during study period.Annual observed rainfall (R) and estimated percolation (P=1.05RÀ ET f À Q À Z r Dh/Dt) are also shown with incomplete years for 2010 and 2016.ET, evapotranspiration.

Fig. 4 .
Fig. 4. (a) Significant power-law relation between monthly deep percolation (P) and streamflow (Q) shown with regression parameters (p = 0.02, q is Spearman's rank correlation coefficient, RMSE is the root mean square error, and n R is the number of months with streamflow).(b) No significant relation between monthly deep percolation (P) and the change in soil moisture storage (Z r Dh/Dt) during the number of months with no streamflow, n NR .

Fig. 5 .
Fig. 5. Comparison of (a) channel and (b) hillslope h at varying depths, along with rainfall (R) and streamflow (Q).No data are available at the channel site prior to 4 August 2016.

Fig. 6 .
Fig. 6.Depth-averaged soil moisture distribution at four time periods during a large storm event on 13 September 2013: (a) antecedent, (b) peak, (c) 48 h after peak, and (d) 120 h after peak.(e) Spatially average soil moisture (h SN ) for event with time periods labeled.

Fig. 7 .
Fig. 7. (a) Rainfall and channel runoff measurements for event on 4 September 2014.(b) Relation between outlet streamflow (Q) and total runoff at two flumes (Flume 2 + Flume 3) normalized by upstream area.(c) Relation between total runoff volume and upstream area over the study period, with inset comparing derived and mapped channel networks.

Fig. 8 .
Fig. 8. Conceptual diagram of deep percolation in channels.Rainfall (R) on hillslopes infiltrates (I) into the shallow soil and overland flow (sheetflow) is generated.Sheetflow interacts with the channel with a storage capacity S c .(a) When sheetflow volume is less than S c (linked to lower R), the channel does not saturate, and deep percolation (P) is observed, but no outlet streamflow (Q) is generated.(b) When sheetflow volume is more than S c , both P and Q are observed.

Fig. 9 .
Fig. 9. Relation between monthly rainfall and percolation (symbols), including a piecewise linear regression (solid lines) with an inflection point at 62.8 mm/ month (dashed line).

Fig. 10 .
Fig. 10.(a) Annual percolation estimates at three rain gauge sites, along with annual rainfall from the Headquarters site.Dashed line represents long-term average P based on Headquarters.(b) Geomorphic regions in the Jornada Experimental Range with locations of rain gauges and the outlets of all first-order watersheds on the piedmont slope derived using a 10m digital elevation model.

Table 2 .
Watershed water balance components and percolation estimates accounting for measurement uncertainty (mm).