Groundwater recharge in desert playas: current rates and future effects of climate change

We conducted both the empirical and modeling components of this study at the Jornada Basin long-term ecological research site (JRN). This site is located near Las Cruces, NM, USA (+32.5 N, −106.8 W, elevation 1188 m) and contains all of the ecosystem types and geomorphic landforms that are typical for systems in the Basin and Range province (Peters and Gibbens 2006). JRN is situated above the Jornada del Muerto aquifer with mean annual precipitation of 237 mm yr⁻¹. (figure 1). The Jornada LTER is composed of two grassland and three shrubland ecosystem types and average temperature of 24 °C. The grassland ecosystems are upland Black Grama grasslands (Bouteloua eriopoda) and lowland playa grasslands co-dominated by Tobosa grass (Pleuraphis mutica), Alkali Sacaton (Sporobolus airoides) and Vine-mesquite grass (Panicum obtusum). The shrubland ecosystems are: Tarbush (Flourensia cernua) on lower piedmont slopes, Creosotebush (Larrea tridentata) on upper piedmont slopes and bajadas, and Honey Mesquite (Prosopis glandulosa) on the sandy basin floor (Peters 2013).

Zoom In Zoom Out Reset image size

There are 100 playas that account for ~1% of JRN area (Gibbens et al 2005). Playa catchments are the upland drainage areas that contain grassland and shrubland ecosystem types. Playas in this region experience periodic flooding (Peters et al 2012) and they are the only known areas that do not contain water-restrictive petrocalic horizons (Havstad et al 2006). Playa soils are well-mixed vertisol soils that are high in clay content and contain Holocene aged lacustrine deposits of the derived mainly from monzonite, rhyolite, and andesite (Wondzell et al 1990).

We used the chloride mass balance (CMB) approach to empirically estimate groundwater recharge rates beneath playa surfaces (Wood and Sanford 1995). CMB is the most common method for estimating groundwater recharge in unsaturated soil zones (Allison and Hughes 1978, Scanlon 1991). Although other more resource intensive tracer methods can be used, CMB is commonly used on its own to estimate groundwater recharge (Scanlon et al 2006, Gates et al 2008, Gurdak and Roe 2010). Scanlon et al (2006) provided many examples of estimating recharge with CMB well above the water table and assuming downward vertical flow of groundwater beneath the rooting zone. This approach was relevant to our study because of the thick (>10 m) unsaturated zone in the Jornada Basin (Havstad et al 2006). Surface inputs of chloride (Cl_p) (mg l⁻¹) and precipitation (P) (mm yr⁻¹) are balanced by the mass out beneath the playa via chloride in the unsaturated zone (Cl_uz) (mg l⁻¹) and deep percolation that results in groundwater recharge (R) (mm yr⁻¹) (equation (1)). We calculated recharge by solving equation (1) for R (equation (1)).

$$\begin{equation} {\rm{R}} = \frac{{\left( {{\rm{C}}{{\rm{l}}_{\rm{p}}}} \right)\left( {\rm{P}} \right)}}{{\left( {{\rm{C}}{{\rm{l}}_{{\rm{uz}}}}} \right)}} \end{equation} \tag{ 1 }$$

We collected one continuous 5 m soil core in the unsaturated zone from the center-low point of each study playa (n = 20). We sampled at the lowest topographic point of each playa, because this area is flooded most-frequently (Gurdak and Roe 2010, McKenna and Sala 2016). We collected each core in August 2014 when playas were not flooded. We took five 100 g soil samples from below the core at 1 m increments (1–5 m). Soils were homogenized and 3 subsamples were taken from each 1 meter depth. We combined gravimetric water content (g g⁻¹) with soil bulk density (g m⁻³) measurements (Elliot et al 1999) to calculate volumetric soil water content. Cl⁻ was measured (mg l⁻¹) using an ion selective electrode. We measured average soil Cl⁻ from across 1–5 m depths for each playa (figure S3). To calculate recharge with the CMB method, we used values of annual precipitation and wet and annual dry deposition values from JRN data. We calculated mean annual precipitation from the 100 year record (1914–2015) of the centrally located Jornada weather station (figure 1). We calculated annual rates of wet and dry Cl_p deposition (mg l⁻¹) using a 30 year record (1983–2013). Wet and dry deposition of Cl_p were measured monthly at the same centrally located weather station using AeroChem Metrics collector (Havstad et al 2006).

Each playa served as an independent unit for our statistical analyses. Under this framework, adding more samples per playa would not have improved the statistical power. Taking several samples per playa and using them as independent replicates would have violated the independence criterion of multiple regression (Kutner 2005).

We used remotely sensed data to measure catchment biophysical characteristics for each of the study catchments. We defined a playa catchment as the upland area that drains into the playa. A playas catchment does not include the playa itself. We used Normalized Difference Vegetation Index (NDVI) data derived from 250 m² resolution satellite Moderate Resolution Imaging Spectroradiometer (MODIS) to calculate average vegetation cover, specifically we used the vegetation continuous fields MOD44B product. Digital soil maps were used to estimate an average soil texture for each playa catchment (Soil Survey Staff 2016). We measured the area and slope of each catchment by analyzing 5 m digital elevation model (DEM) data from http://jornada.nmsu.edu/lter/data/spatial. We used multiple regression analysis to assess correlations between biophysical catchment characteristics and groundwater recharge in playas. All analyses were conducted using R version 3.0.2 (R 2016). The best-fit-model was chosen using Akaike information criterion (AIC) (Kutner 2005).

We used the Limburg Soil Erosion Model (LISEM) (De Roo et al 1996) to simulate playa runon for 20 playa catchments over a 20 year period. LISEM is a well-known hydrologic model that has been successfully used to simulate runoff for a variety of catchments around the world (Cuomo et al 2015). This model has been used to simulate runoff during and immediately after rainfall events in dryland catchments from 0.1–100 km (De Roo and Jetten 1999, Hessel et al 2006, Baartman et al 2012). Runoff, infiltration, and interception were calculated in LISEM using data from spatially distributed soil, vegetation, and elevation maps for each catchment (De Roo et al 1996). Runoff occurs as infiltration excess Hortonian overland flow (Baartman et al 2012). Dry streambeds act as preferential flow paths for delivering upland runoff to playas (Havstad et al 2006). Catchment maps were created using PCRaster (Schmitz et al 2014). We calculated soil-physical characteristic using soil texture from (Soil Survey Staff 2016), and literature values (Rawls et al 1983). Saturated conductivity (mm hr⁻¹), soil water tension (cm), and saturated volumetric soil moisture content for each soil type were derived from literature (Rawls et al 1989). Initial volumetric soil moisture values were gathered from monthly JRN soil moisture values (1979−1989) (Nash et al 1991). Literature values were used for Manning's surface resistance (n), random roughness coefficients, and vegetation height (cm) (Weltz et al 1992). We modeled playa runon using JRN hourly precipitation record 1992–2011. Throughout those 20 years, there were 560 unique rainfall events above 1 mm in size (figure S1 available at stacks.iop.org/ERL/13/014025/mmedia). To validate the model, we compared outputs of runon to observed playa flood volume for 14 rainfall events (figure S2). We then ran the model for all 560 rainfall events from 1992–2011.

Zoom In Zoom Out Reset image size

After modeling runon from recent-past precipitation events, we independently evaluated the effect of increased precipitation variability and decreased mean annual precipitation. We used predictions of both precipitation variability and amount from three different CMIP5 representative concentration pathways (RCPs) (Wuebbles et al 2014). We manipulated the 560 rainfall events from our historical 20 year record and generated a new 560 event series for RCP 4.5, RCP 6.0, and RCP 8.5 scenarios. We chose a range of emission scenarios to document the largest potential effects on rainfall-runoff-recharge rates. Using modeled relationships from Sun et al (2007), we decreased 0–10 mm events (n = 406) by 2% per °C increase. We increased 10–20 mm events (n = 99) by 5%, 20–50 mm events (n = 53) by 6%, and >50 mm events (n = 1) by 7% per °C increase. These changes maintained mean annual precipitation constant and increased the standard deviation and mean event size (table S1). For this region, global circulation models predict precipitation amount to decrease 2% per °C increase in atmospheric temperature (Pierce et al 2013). To simulate a decrease in mean precipitation for RCPs, we reduced the size of all precipitation events accordingly. We used regression analysis to determine the relationships between playa groundwater recharge (mm yr⁻¹) and modeled runon (m yr⁻³), as well as the playa runon (m event⁻³) and precipitation event size (mm event⁻¹). We used the rainfall-runon and runon-recharge regression models to calculate how changes in both precipitation variability and mean precipitation amount affect recharge underneath playas.

Our empirical estimates of recharge showed evidence of groundwater recharge in the unsaturated zone beneath playas (figure 1) (table S2 and figure S2). We estimated groundwater recharge occurring beneath 100% of the sampled playas (figure 1). The average recharge of playas was 6 mm yr⁻¹, which was 2% of the annual rainfall for the study area. Recharge beneath playas ranged from 0.10–28 mm yr⁻¹, and 65% of the playas sampled had recharge rates below 5 mm yr⁻¹ (figure 1). Recharge rates higher than 15 mm yr⁻¹ occurred only in 15% of playas (figure 1).

We found that differences in groundwater recharge rates among playas were correlated to the size, slope, and soil texture of each catchment (R² = 0.78, p < 0.001, AIC = 53.6). Recharge rate increased with area and slope, and decreased with percent sand of a catchment (figure 2). We found the best-fit model for predicting groundwater recharge to be: Groundwater recharge (mm yr⁻¹) = 29.82 + 0.27⁎catchment area (ha) + 0.82⁎Catchment slope (% rise) −0.38⁎catchment soil texture (% sand). Catchment vegetation cover was not significantly correlated to playa groundwater recharge. These results suggest either vegetation cover did not physically affect runoff production or the range of vegetation cover in catchments existing in our study site was not large enough to capture the physical effect of vegetation on runoff production.

We found that there was a highly significant (R² = 0.77, p < 0.05) linear relationship between the simulated amount of annual runon a playa received and the observed amount of annual groundwater recharge that occurred beneath that playa: Groundwater recharge = −0.23+ 0.0021⁎runon (m yr⁻³) (figure 3). From this relationship and our earlier correlations between catchment area, slope, and soil texture, we can also infer that a playa with a combination of the largest area, steepest slope, and least sand (figure 2) would produce the most runon for the adjacent playa. In order to estimate how changes in precipitation would affect groundwater recharge, we analyzed how the size of individual precipitation events controlled playa runon. Playa runon occurred when precipitation events were larger than 20 mm, which happened twice yearly on average, and only 8% of all rainfall events over a 20 year period were above 20 mm (figure 4). When rainfall events were above 20 mm, the size of a rainfall event correlated linearly (R² = 0.81, p < 0.001) to the average amount of playa-runon generation: Mean playa runon = −6369.86 + 259.78⁎rainfall event size (mm) (figure 4). We used our empirically derived rainfall-runon (figure 4) and runon-recharge (figure 3) relationships to assess how projected changes in precipitation would influence runon and subsequently change groundwater recharge rates in playas.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

After modeling runon from recent-past precipitation events, we independently evaluated the effects of increased inter-annual precipitation variability and decreased annual precipitation amount on playa groundwater recharge. We manipulated a 20 year rainfall record to reflect both increased precipitation variability and decreased precipitation amount. We then used our rainfall-runon-recharge models (figure 3 and figure 4) to calculate the response of groundwater recharge beneath playas to changes in precipitation variability and amount under RCP 4.5, RCP 6.0, and RCP 8.5 scenarios (figure 5). Climate models predict that increased atmospheric temperatures will increase precipitation variability and decrease precipitation amount in Southwestern USA (Melillo et al 2014). Under different climate-change scenarios for Southwestern USA, average atmospheric temperature would increase between 2 and 6 °C for RCP 4.5 and RCP 8.5 (table S1). A rise in temperature would increase precipitation variability between 5%–17% depending on the scenario (Wuebbles et al 2014). In turn, increased variability will increase the frequency of occurrence of large (>20 mm) rainfall events and the mean event size (Sun et al 2007) (table S1). We found that for every 1% increase in inter-annual precipitation variability, average playa groundwater recharge rates increased 18%. In the most-extreme scenario, average playa groundwater recharge rates increased 300% from 6 mm yr⁻¹ to 22 mm yr⁻¹ (figure 5).

Zoom In Zoom Out Reset image size

Climate change predictions for Southwestern USA call for a decrease in mean annual precipitation of 2% for every degree Celsius increase in temperature (Wuebbles et al 2014). Mean annual precipitation would decrease 4%–12% under the RCP 4.5–RCP 8.5 scenarios (Pierce et al 2013). We found that for every 1% decrease in precipitation amount, average playa groundwater recharge rates decreased 5%. In the extreme scenario, average playa groundwater recharge rates decreased 50% from 6 mm yr⁻¹ to 3 mm yr⁻¹ (figure 5). Overall, we found that climate change would have a net positive effect on playa groundwater recharge resulting mainly from an increased number of large runoff-generating rainfall events.