Climate and topographic controls on simulated pasture production in a semiarid Mediterranean watershed with scattered tree cover

Natural grasses in semiarid rangelands constitute an effective protection against soil erosion and degradation, are a source of natural food for livestock and play a critical role in the hydrologic cycle by contributing to the uptake and transpiration of water. However, natural pastures are threatened by land abandonment and the consequent encroachment of shrubs and trees as well as by changing climatic conditions. In spite of their ecological and economic importance, the spatiotemporal variations of pasture production at the decadal–century scales over whole watersheds are poorly known. We used a physically based, spatially distributed ecohydrologic model applied to a 99.5 ha semiarid watershed in western Spain to investigate the sensitivity of pasture production to climate variability. The ecohydrologic model was run using a 300-year-long synthetic daily climate data set generated using a stochastic weather generator. The data set reproduced the range of climatic variations observed under the current climate. Results indicated that variation of pasture production largely depended on factors that also determined the availability of soil moisture such as the temporal distribution of precipitation, topography, and tree canopy cover. The latter is negatively related with production, reflecting the importance of rainfall and light interception, as well as water consumption by trees. Valley bottoms and flat areas in the lower parts of the catchment are characterized by higher pasture production but more interannual variability. A quantitative assessment of the quality of the simulations showed that ecohydrologic models are a valuable tool to investigate long-term (century scale) water and energy fluxes, as well as vegetation dynamics, in semiarid rangelands.


Introduction
Traditional Mediterranean agrosilvopastoral systems support high levels of biodiversity in a wide variety of coexisting natural and man-made habitats, such as grazing areas, agricultural lands, scrublands, forests or wildlife spaces Campos-Palacín, 2004). Natural grasses and pastures are an important element of cohesion between these habitats by supporting livestock and other fauna, by protecting the soil against erosion and degradation, and by controlling the soil hydrologic and thermal regime (Schnabel, 1997;Paço et al., 2009). The economic importance of pasture encourages the proper management and conservation of Mediterranean agrosilvopastoral systems; however, owing to climate characteristics of semiarid Mediterranean environments, natural herbaceous production is highly variable with a pronounced seasonality, being highest in spring, low in autumn and winter, and nil during summer (Montero et al., 1998;Joffre and Rambal, 1993). Additionally, pasture yield is usually low and its spatiotemporal distribution is strongly conditioned by the balance of positive and negative effects of limiting factors such as water, light, or nutrients (Brooker et al., 2008).
Decreased pasture yields may upset the balance of habitats and threaten the sustainability of these Mediterranean systems due to changes in land use associated with a revision of economic priorities and management decisions. Indeed, pastures in Mediterranean Europe have been experiencing land abandonment and consequent encroachment of shrubs and forest (Rivest et al., 2011;García-Ruiz and Lana-Renault, 2011;Lavado-Contador et al., 2004), which may lead to increased competition for resources, such as water and light, among different layers of vegetation (Cubera and Moreno, 2007a). The abandonment of traditional agrosilvopastoral systems may not only have important ecologic consequences but may also have a significant impact on regional economies and on food security by affecting forage quality and quantity and by affecting productivity and protection of the agricultural landscape against degradation.
Improved knowledge of the frequency of low and high pasture productivity periods and the expected variability of yields in different locations of a region permits making better informed management decisions that contribute to the sustainability of agrosilvopastoral systems; however, we still only have a partial understanding of the ecohydrological processes that control plant productivity across space and time (Asbjornsen et al., 2011).
From the mid 90's there has been a growing interest in the complex interactions between ecological and hydrological processes at multiple scales (Viville and Littlewood, 1996;Rodríguez-Iturbe, 2000;Wang et al., 2012;Caylor et al., 2005;Caylor et al., 2009;Porporato et al., 2002;Rodriguez-Iturbe et al., 1999). Because of the complex and non-linear interactions between vegetation and hydrology, few studies focus on the larger scales, such as landscapes or watersheds, where the processes are less understood (Asbjornsen et al., 2011). A limited number of models have been developed in the last decade to investigate ecohydrologic interactions at watershed and regional scales (e.g., Ivanov et al., 2008;Oleson et al., 2010;Tague and Band, 2004;Maneta and Silverman, 2013;Fatichi et al., 2012). Most of the studies using these models have focused on short-term studies because of the long run times derived from their complexity and because the lack of existing extensive climate data sets (longer than a few decades) needed to force the models. These limitations have resulted in few studies conducting simulations over the entire range of ecohydrological conditions that can be expected under current climate variability. These studies would be highly valuable to improve our understanding of the variability of pasture production and to inform grassland management.
Reproducing the entire range of ecohydrologic states to capture relevant watershed processes requires the ability to simulate extensive periods in the order of hundreds of years at small spatial (1-50 m) and temporal (daily) scales. Maneta and Silverman (2013) present a ecohydrologic model with a level of complexity that can make the simulation of extensive periods at detailed spatial and temporal scales tractable while maintaining a strong mechanistic description of the processes. The lack of extensive input data sets to the model can be overcome by producing synthetic data sets with stochastic weather generators (SWG). These tools have been successfully used since the early 1980s (Richardson, 1981) to generate long time series of synthetic weather data that are statistically indistinguishable from observed shorter term climate records (Semenov and Barrow, 2002). SWGs have been used to simulate future scenarios of climate change (Fatichi et al.,   2011; Semenov and Barrow, 1997), crop yields (Semenov and Porter, 1995;Ivanov et al., 2007) or regional hydrologic response (Xia, 1996;Dubrovský et al., 2004).
In this paper we use a combination of mechanistic models and SWG to investigate the spatiotemporal variability of pasture production at watershed scales relevant for management. Questions that we seek to address include: how does pasture production respond to climate variability in combination with antecedent basin conditions? How sensitive is the production of pasture to the temporal distribution of precipitation during the year? How important are topographic controls vs climatic controls in determining the spatial and temporal dynamics of production in a watershed? Does the relative importance of these controls vary for different years and under different circumstances?
While abundant studies have applied numerical models to the study of grassland productivity (Montaldo et al., 2005;Istanbulluoglu et al., 2012) and some work has a focus on the spatiotemporal variability of pasture production over long periods (century scale) and large areas (Clark et al., 2003;Tubiello et al., 2007), to the authors' knowledge no studies have applied comprehensive mechanistic numerical models to address the questions posed above. Experimental or field studies have not addressed satisfactorily these questions either because pasture production over large areas is typically determined with a limited number of measurements commonly taken over a few years and at very specific locations (Plaixats et al., 2004;Santamaría et al., 2009). The limited number of samples could provide a skewed or erroneous estimate of the actual long-term pasture production of a region or farm because short-term studies with infrequent sampling may not properly capture the effect of weather variations, such as wet and dry periods, and the specific sampling locations may not properly characterize the actual spatial variation. A modeling approach is therefore preferred in this study.

General description
The study area is an experimental drainage basin located in the southwestern part of the Iberian Peninsula with an area of 99.5 ha (Fig. 1), characterized by an agrosilvopastoral land use system called dehesa in Spain. Geologically, the study area forms part of the Iberian Massif of Precambrian age, where the dominant rocks are greywacke and schist, which were eroded giving rise to an erosion surface. Topography of the drainage basin is gently undulating with an average elevation of 394 m a.s.l., where SSW is the dominant aspect. The climate is Mediterranean with a high seasonal and interannual rainfall variability , which determines the available water content for plants, and a marked dry season during summer that can last four months or even more. Average annual precipitation for the period between 1999 and 2012 was 488 ± 149.5 mm (mean ± standard deviation) and mean monthly temperatures ranged between 7.4 ± 1.7 • C in January to 26.4 ± 1.5 • C in July and August. Annual potential evapotranspiration is twice the annual rainfall amount. Vegetation is typically Mediterranean, characterized by a two-layered vegetation structure, with a layer of scattered trees (Quercus ilex) at low density (20 ± 18 individuals ha −1 ), and a pasture layer. Natural pastures are composed of annual and perennial herbaceous plants, abounding especially annual grasses (such as Vulpia bromoides, Bromus sp. or Aira caryophyllea) and annual legumes (Ornithopus compressus, Lathyrus angulatus and several species of Trifolium), starting to grow with the first rainfall in autumn and reaching maximum production in spring. A layer of shrubs is also frequent (Retama sphaerocarpa), commonly eliminated by ranchers to facilitate pasture growth.
Soils in the catchment have a high bulk density (≈ 1.5 g cm −3 ) are poor in nutrients and have low organic matter conten (≈ 3 %) except below tree cover where it is higher in the upper 5 cm (Schnabel et al., 2013b). Roots are concentrated in the upper soil layer (Moreno et al., 2005), favoring the higher porosity (≈ 45 %) of the topsoil. Two geomorphologic units can be distinguished in the catchment which determines the type of soil and its hydrologic properties. The boundary between these units is marked by the 395 m contour (Fig. 1). The geomorphological unit above 395 m is the northern part of the catchment. It constitutes the slopes of a pediment with sandy loam soils classified as Luvisols (FAO, 1988), rich in rock fragments that provides it with a higher permeability and saturated hydraulic conductivity than the remaining soils (Van Schaik et al., 2008;Van Schaik, 2009). Soil depths in this unit are variable, often exceeding 1 m to bedrock and with an argillic B horizon. The other geomorphologic unit, flat to gently undulating, is located in the lower part of the basin. In this unit soils are very shallow (Cambisols and Leptosols), ranging between 20-50 cm, developed on impervious bedrock of schist and greywacke, which frequently outcrops. The lowest areas of this unit correspond with valley bottoms covered by alluvial sediments reaching a thickness of approximately 1 m in areas next to channels. The main channel is incised into these sediments, actively eroding at present and can be classified as a gully (Gómez-Gutiérrez et al., 2009). Owing to low permeability of these layers some sites are prone to ponding in wet periods (Cerdá et al., 1998;Van Schaik, 2009), which provide an extra water storage that may lengthen the phenological period of the herbaceous plants and that is totally dried in summer. A complete and detailed description of the study area can be found in Maneta (2006) and Van Schaik (2010).

Meteorological data
The study area is equipped with a meteorological station that collects information on precipitation, temperature, relative humidity, global radiation, net radiation, wind speed and direction at intervals of 5 min since the year 2000. Rainfall is also measured in five other locations ( Fig. 1) with tipping bucket type rain gauges of 0.2 mm resolution. This information was aggregated in daily intervals for this study.

Soil moisture content and soil temperature
Volumetric soil water content was monitored by capacitive sensors (Decagon Device, Inc. model EC-5) at 5, 10, 15 and 30 cm depth every 30 min. Soil temperature was measured at 5 cm depth near the soil moisture probes (Decagon Device, Inc. model RT-1). The accuracy of the soil moisture sensors was improved by calibration following the method of Cobos and Chambers (2010). The sensors were grouped in soil moisture stations (SMS) at two sites: site 1 representative of hillslopes with Luvisols, and site 2 representative of the lower part of the catchment with shallow soils. A third SMS was installed in the eastern part of the catchment (Fig. 1). The selection of sites to install the SMSs were based on previous studies by Lavado-Contador et al. (2006), Maneta et al. (2007Maneta et al. ( , 2008a and Van Schaik et al. (2008). The SMSs in site 1 and site 2 began to register in March 2009, while SMS-3 started in May 2010. In each site there are sensors in open grass areas and under tree canopies. The overall soil moisture of each site was considered to be the depthaveraged soil moisture of the sensors under trees and in open areas, weighted by the relative canopy cover in its pixel.

Pasture production
We have measured natural pasture production at site 1 and site 2 for three hydrologic years (from September 2008 through August 2011). To prevent grazing, twelve 1 m × 1 m livestock exclusion cages were installed at midslope positions in open space. Only aerial (above-ground) production is considered in this study. Grasses and forbs were cut twice a year (at the end of winter and at the end of spring), dried during 48 h in an oven at 105 • C and weighed to determine aerial dry matter (DM) production (kg DMvha −1 ).
Measurements of DM were augmented with measurements of pasture height. At each SMS, 16 measurements of plant height were taken biweekly during two hydrological years (from 1 March 2011 to 31 August 2012). The pasture production database was extended by estimating DM from pasture height measurements using their allometric relationship (r 2 = 0.68, n = 12).

Ecohydrologic model
To simulate water and energy exchanges and pasture production, we used a spatially distributed ecohydrologic model as described in Maneta and Silverman (2013). This model couples a two layer (canopy and understory) vertical local closure energy balance scheme, a hydrologic model and a carbon uptake and vegetation growth component. The model was run using climate information from a stochastic weather generator as described below.
Vertical energy transfers are calculated using first-order closure profile equations for momentum, heat and mass under neutral stratification based on flux gradient similarity (Arya, 2001;Foken, 2008). The energy balance is solved for the canopy layer and then for the soil layer using canopy temperature and soil temperature as the closure variables, respectively. Canopy conductance is calculated with a Jarvistype multiplicative model (Cox et al., 1998;Jarvis, 1976). The model takes into account the vertical and lateral redistribution of water and considers the effect of topography. Water can infiltrate into the soil or become runoff, which can reach the channel and exit the watershed, or re-infiltrate downslope. Water infiltration is calculated using the Green and Ampt approximation to Richard's equation (Chow et al., 1988). Lateral water transfers in the soil are simulated using a 1D kinematic wave model (Singh, 1997). Infiltration and lateral subsurface flows are controlled by soil hydraulic properties (hydraulic conductivity, porosity) and by the topographic gradient. The bedrock at the bottom of the soil is considered to be impermeable and when the soil is fully saturated, return flow occurs. Interception of water by canopies is simulated using a bucket model. The forest growth and carbon uptake components are based on 3-PG (Landsberg and Waring, 1997); see Maneta and Silverman (2013) for further details.
The ecohydrologic model by Maneta and Silverman (2013) was extended in this study with a new grass growth component. Net primary production of grass is related to the available radiation intercepted by the canopy and the water transpired: where NPP is net primary production, PAR is photosynthetically active radiation intercepted by the canopy, Transp is transpiration, α is a constant light use efficiency parameter, β is a constant water use efficiency parameter, f (T a ) is a production efficiency function dependent on air temperature (Landsberg and Waring, 1997), and C NPP is a GPP to NPP conversion factor. Transpiration is calculated from the latent heat term of the energy balance equation for the canopy layer, which takes into account relevant environmental conditions (e.g., air temperature, vapor pressure deficit, soil moisture). Aerodynamic resistance and interception of PAR are related to the leaf area index of vegetation as described in Maneta and Silverman (2013). The onset of the growing season and the initiation of dormancy are determined by a threshold in the minimum daily air temperature. NPP is allocated to two carbon pools: aboveground biomass (leaves) and belowground biomass (roots). Aboveground biomass is further divided into green aboveground biomass and dead aboveground biomass. The dynamics of these carbon pools are described by three ordinary differential equations that track their mass balance (Montaldo et al., 2005;Istanbulluoglu et al., 2012): where M g , M r and M d are dry mass in the green grass, root, and dead grass pools, respectively; k sg , k sr and k sd are constant decay coefficients for green, root and dead biomass, respectively. Parameter ξ sd is an adjustment factor for the coefficient of dead biomass decay. This adjustment permits to account for reduced decay during the cold season when the temperature of the canopy (T c ) drops below a given temperature threshold (T ξ ): Parameter a (Eq. 2a, b) controls the allocation of NPP to the aboveground (green leaves) and belowground (roots) pool of carbon based on the spare capacity of the land to carry aboveground biomass (Istanbulluoglu et al., 2012) : where LAI g , LAI max , and LAI d are green, maximum, and dead grass leaf area indices, respectively. The denominator of Eq. (4) indicates the space available to grow green leaves.
The transformation of the aboveground mass to leaf area index is done using the specific leaf area index for green and dead leaves: where σ LAI g and σ LAI d are the specific leaf area indices for green and dead leaves. Total leaf area index (LAI t ) is considered to be the sum of the green and dead leaf area indices.

Model setup Hydrologic properties, land cover and vegetation parameters
The modeling domain was discretized with a 30 m × 30 m grid, as used in previous studies (Maneta et al., 2008). A digital elevation model (DEM) was used to delineate the limits of the basin, obtain a map of local slopes and other basic information on the geometry of the domain. The drainage direction network was calculated using a deterministic steepest descent algorithm (D8 algorithm). Maps of soil properties such as soil depth, porosity, and other hydrologic properties ( Fig. 2) where derived from the geomorphologic characteristics of the basin as described in Maneta et al. (2008). Soil albedo, emissivity and soil thermal capacity were considered uniform in space. Tree density and tree canopy cover maps were obtained manually digitizing a point for each individual tree in a high-resolution aerial photograph, then calculating the density of points using a 3 × 3 moving average kernel. The fraction of the area covered by canopy was calculated using a maximum likelihood supervised classification technique from a 24 bit color submetric-resolution aerial photography. Once a canopy mask was produced, the canopy coverage was obtained by calculating the fraction of pixel classified in each of the larger pixels used in the simulation ( Fig. 2) (Maneta, 2006). Physiological and structural parameters for trees (Quercus ilex) were taken from the literature (Table 1), while parameters related to pasture were mostly manually adjusted (Sect. 3.4).

Generation of atmospheric forcing
LARS-WG v5.5 (Semenov and Barrow, 2002) is a SWG that generates temporal series of synthetic weather statistically similar to observations at a single site. LARS-WG generates the synthetic weather by sampling from semi-empirical distributions that takes into account the length and the frequencies of wet and dry periods and the covariance among variables, which is important to properly simulate Mediterranean climates. More information about this SWG can be found in Semenov et al. (1998). accumulation (number of pixels that spill on another), E) tree density (trees ha -1 ), F) tree canopy cover (0-1). Maps were obtained as described in Maneta et al. (2008). We used 13 years of data from our meteorological station (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012) to inform LARS-WG about weather patterns in our basin. We assume that the 13 years of available data are representative of the current climate. Small gaps in the data set were filled using data from a meteorological station located at a distance of 24 km from the study area. A linear regression model relating data between the stations was sufficient to correct satisfactorily the differences in the external station. LARS-WG was applied to generate a series of 300 years of minimum and maximum temperature, precipitation and solar radiation at the daily timescale. The generation of a 300-year-long climate data set was chosen to ensure that we are capturing the most common combinations of weather events and basin antecedent conditions that ranchers are likely to experience during the growing season. Other atmospheric information necessary to run the model was generated as follows: daily relative humidity was estimated with a multiple regression model that used daily mean, maximum and minimum temperature and daily rainfall as predictors (r 2 = 0.75). Wind velocity was obtained by repeating a series of 51 years extracted from a station located  Infante et al. (2003) at 24 km from the study site. Daily long wave radiation was estimated from air temperature using the method described by Swinbank (1964).

Model calibration, spin up and data analysis
The calibration runs were done running the period from 1 September 2008 to 31 August 2012 in a continuous loop using daily time steps. Model parameters listed in Table 2 were manually calibrated until soil moisture, soil temperature and pasture yield achieved steady state and satisfactorily matched the available measurements of soil moisture, soil temperature, and pasture yield based on height measurements. Calibration was based on trial and error systematically changing parameters one at a time. When available, the initial trial value was based on values cited in the literature or based on experience. Model performance was quantified using the coefficient of determination, root mean square error, bias and Nash-Sutcliffe efficiency coefficient between modeled and observed soil moisture, soil temperature and pasture yield.
Once performance was satisfactory with parameter values within a realistic range the model was considered calibrated. The calibrated model was used in a 300-year-long simulation at daily time steps resulting in 109 500 maps per state variable reported by the model. State variables analyzed included soil moisture, soil temperature, pasture production, pasture evaporation and transpiration, and tree evaporation and transpiration. Time averages and standard deviations for the entire simulation period were calculated for each variable, except for pasture production. For this latter variable, the average and standard deviations for 1 June were used in the analysis because this date corresponds to the end of the vegetative period of herbaceous plants and can be considered as the day of maximum accumulated production.

Model performance
Mean annual precipitation for the simulated period was 508.8 mm with a standard deviation of 118.2 mm. Maximum and minimum annual rainfall were 934.1 mm and 188.2 mm, respectively. The longest dry spell spanned four years with annual rainfalls lower than 386.9 mm year −1 , while the maximum wet period lasted three years with rainfall in excess of 693.4 mm year −1 . A comparison between simulated and observed atmospheric data indicated that the SWG was properly calibrated and that it successfully generated a synthetic times series that was statistically indistinguishable from the observations (Table 3) except for rainfall in July and August. This is because during these months precipitation volumes are insignificant and small fluctuations about the very low observed precipitation values have a relatively large influence in the K-S statistic. This is of minor importance because rainfall in these months is virtually zero. Further inspection of the results showed that the generated weather series represents the seasonal and interannual variations typical of the Mediterranean climate.
An initial inspection of the graphs shown in Figs. 3 and 4 indicates that the model reproduced (to a high degree) the  observed dynamic of soil moisture and temperature. The simulation captured the seasonal variations of soil moisture, including the wetting and recession rates, but also much of the observed high-frequency variation. Some mismatch can be observed in the reproduction of wetting peaks, such as   those of site 1 (Fig. 3a). There is a general dampening of the amplitude of high-frequency variations that may be due to the model representation of soil moisture as the average over the entire soil profile (Maneta and Silverman, 2013). However the standard goodness-of-fit statistics and descriptive statistic confirmed a satisfactory fit with high coefficients of determination (r 2 ≥ 0.80), low RMSE (≤ 0.047 m 3 m −3 ) and similar statistics for all measurement stations (Table 4). Further evaluation of the model performance show high Nash-Sutcliffe coefficients (≥ 0.75) and low prediction bias (≤ 0.018 m 3 m −3 ).
The simulated soil temperature captured the highfrequency variation of observed soil temperature (Fig. 4). However, during the first year simulated temperatures were higher than observed in both study sites, which could be caused by uncommonly low pasture yields simulated that year and hence an overestimation of the amount of radiation reaching the bare soil, while actual ground covered by pasture was much higher at the SMS sites because they were protected against grazing. Efficiency statistics for soil temperature were satisfactory, with coefficients of determination r 2 ≥ 0.89 and the Nash-Sutcliffe efficiency criterion 0.86, increasing our confidence on the capacity of the model to represent the energy fluxes in the study site (Table 4).
Simulated annual pasture production matched well the observed data at both field sites ( Table 4). The average simulated value of production for both sites was 630.9 kg DM ha −1 , very similar to the observed 623.8 kg DM ha −1 . Other descriptive statistics (minimum, maximum, standard deviation) and goodness-of-fit statistics confirming the model in our research area are shown in Table 4. The model produced a satisfactory description of the spatiotemporal dynamics of production, which is supported by the high prediction efficiency of the model (Nash-Sutcliffe ≥ 0.75; r 2 ≥ 0.76) and low residual errors (RMSE = 164.8 kg DM ha −1 ).
The phenological cycle of the herbaceous plants in the study site (Fig. 5) is captured in the simulated data and includes low production in autumn although dependent on antecedent precipitation, scarce production in winter because of low air temperatures and available energy, high production in spring when water and energy are available and an absence of production in summer because of lack of water. It is important to note that once pasture is cut at the sites to measure its dry biomass, the exclusion cage is moved to a nearby location, which contributes to the difference between DM estimated from cuts (blue diamonds) and from vegetation height (green circles) since production is highly variable even at short distances (as indicated by the standard deviation of pasture cuts, Fig. 5). In contrast, plant height is always and consistently measured at the same location (SMS). Even though we do not have direct measurements of tree transpiration to verify our simulations, it is of value to compare our results with the transpiration of Q. ilex reported in the literature. Figure 6 shows tree and pasture transpiration during four hydrological years in a pixel of site 1 and site 2. Simulated dynamics of tree transpiration in site 1 follow a marked seasonal cycle reaching maximum values in spring when environmental conditions were optimal for growth. The maximum simulated value was 1.0 mm d −1 which is slightly lower than observed values reported by Infante et al. (2003), who measured maximum daily transpiration between 1.2 and 1.4 mm d −1 . Higher values were found by Paço et al. (2009), who even observed values exceeding 2.5 mm d −1 . Q. ilex maintained transpiration throughout the year, even during summer when the soils are dry.
Pasture transpiration is associated with the seasonal phenological cycle typical of annual herbaceous plants. In both sites, low transpiration occurred in autumn and is associated with low pasture growth (Fig. 6). Maximum values were registered in spring, not exceeding 1.75 mm d −1 , when herbaceous plants find the most suitable environmental growth conditions. Similar values were also observed by Paço et al. (2009) in an analogous ecosystem, where the authors estimated maximum peaks in excess of 1.5 mm d −1 , while Joffre and Rambal (1993) found different values depending on the annual rainfall in more humid dehesas, ranging from 2.0 to 2.9 mm d −1 .

Spatial distribution of soil moisture and evapotranspiration
Simulated average catchment soil moisture for the 300 years was 0.158 m 3 m −3 , although strong variations were found among different locations in the study area ranging from 0.070 to 0.285 m 3 m −3 (Fig. 7a). Average simulated soil moisture at site 1 was slightly lower than at site 2, with 0.174 and 0.201 m 3 m −3 , respectively, which is in accordance to the observed differences between sites of measured values (Table 4). A multiple regression analysis revealed that the most explanatory variables determining the spatial distribution of soil moisture are canopy cover, porosity, slope, and elevation. These variables explained 68 % of the observed variance and, with the exception of porosity, showed a negative correlation with soil moisture. Canopy cover showed a particularly strong negative relationship with soil moisture, indicating that the reduction of water reaching the ground due to rainfall interception and the additional water uptake by the trees was a more determinant control of soil moisture than the reduction of incident radiation and evaporation below tree canopies due to shading.
Low lying areas had greater average soil moisture (Fig. 7a). These areas correspond to the valley bottoms and flat footslopes, which show better conditions for water maintenance by the effect of topography (concentrating water) or thicker soils with a higher content of clay and silt particles and greater porosity (McGlynn et al., 2003;Jencso et al., 2009). In contrast, hillslopes and areas at greater altitude had lower soil moisture values, which could be attributed to smaller contributing areas, higher canopy cover and coarser J. Lozano-Parra et al.: Climate and topographic controls represents simulated average pasture yield for whole pixels in every Site, with +/-1 standard deviation (green shade). Green circles represent average pasture production based on height measurements; blue rhombuses represent average pasture production based on plant cuts (moustaches correspond to +/standard deviation).

Fig. 5.
Observed and simulated accumulated pasture production at (A) site 1; and (B)site 2. The red line represents simulated average pasture yield for whole pixels in every site, the green shade represents the ±1 standard deviation of spatial variability computed from the 8 neighboring cells. Green circles represent average pasture production based on height measurements; blue diamonds represent average pasture production based on plant cuts (error bars correspond to ±1 standard deviation of local spatial variability of pasture). soil textures. However, a small area in the northeastern upper part of the catchment also showed high average soil moisture values, which could be explained by its low tree density and low canopy cover.
These results highlight the importance of trees in the spatial distribution of soil moisture. This has been observed in dehesa systems by Lavado-Contador et al. (2006), Martínez Fernández et al. (2007 or Moreno and Cubera (2008). Whether trees enhance or reduce soil moisture with respect to open areas seems to be dependent on the climatic conditions of the site (Lozano-Parra et al., 2011). Joffre and Rambal (1988) found higher water content beneath tree canopies in sub-humid ecosystems, which could explain enhanced pasture yields in these situations. Likewise, Gindel (1964) observed also higher water content beneath canopy than in open areas under subtropical and semi-desert conditions. In contrast, García-Estringana et al. (2013) measured lower soil  moisture under forest cover in a Mediterranean mountain area, while Cubera and Moreno (2007b) and Gea-Izquierdo et al. (2009) found lower water contents beneath canopy in semiarid conditions with scattered trees, which is in accordance with our results.
The variability of soil moisture is presented in Fig. 7b and shows a spatial distribution that correlates with the distribution of soil moisture averages. Higher temporal variability of soil moisture was observed in areas with high average soil moisture (e.g., valley bottoms). In contrast, areas with low mean soil water content such as hillslopes with high gradients showed less temporal moisture variability. An explanation for this behavior is that regions with intermediate and higher water contents and soils with good retention properties have more opportunities for soil moisture fluctuations than drier soils with poorer soil water retention capabilities that quickly drain and dry.
Simulated evapotranspiration was marked by the spatial distribution of vegetation cover and by topography (Fig. 7c). Maximum values were found in the valley bottoms where water content remained high during most of the year. High values were also observed in areas with high tree density, while they were lower in open areas where herbaceous vegetation dominates. Annual mean value of actual evapotranspiration for the whole catchment was 390 mm while annual mean precipitation was 508 mm. This implies that about 120 mm could become runoff or to be stored in the soil reservoirs ( Fig. 1) or rock fractures of the impermeable bedrock of the catchment. In support of this, Schnabel et al. (2013a) measured in the same environment runoff values that oscillated between 10 and 190 mm depending on annual precipitation. The simulated annual evapotranspiration values in areas of relatively high tree density are similar to the 590 mm reported by Joffre and Rambal (1993) under tree cover in sub-humid Mediterranean rangelands. They found, however,  higher annual values, 400 mm, in open spaces, which could be explained because their study was carried out in a wetter environment.

Pasture production: temporal dynamics
At site 1 annual average dry matter production was 338.0 kg ha −1 , with a standard deviation of 172.5 kg ha −1 , and maximum and minimum values of 977.6 and 20.7 kg ha −1 year −1 , respectively (Table 5). At site 2 annual average dry matter production was higher (456.0 kg ha −1 ), also with higher maximum (1030.9 kg ha −1 year −1 ) and minimum (29.9 kg ha −1 year −1 ) values of annual dry matter production. Site 1 showed higher relative variation of production Figure 8. Simulated average pasture production and precipitation at Site 1, A) at the annual timescale for 300 years. B) for ten years at the daily timescale (the green shade represents +/-1 standard deviation of pasture production, and the blue bars is the rainfall) Fig. 8. Simulated average pasture production and precipitation at site 1; (A) at the annual timescale for 300 years; (B) for ten years at the daily timescale (the green shade represents ±1 standard deviation of pasture production, and the blue bars is the rainfall). as compared to site 2. Coefficients of variation for each site were 0.51 and 0.40, respectively.
Also, the range of pasture production was slightly higher at site 2 (approximately 1000 kg DM ha −1 year −1 compared to 957 kg DM ha −1 year −1 for site 1). These production values rank the study site as a low productivity rangeland that requires the introduction of supplementary fodder to maintain livestock. Bell (2006) reports that the critical pasture mass necessary to sustain a sheep ranch is between 400 and 1700 kg DM ha −1 , while for cattle 700 to 2900 kg DM ha −1 . Productivity values for similar Mediterranean rangelands are highly variable, as reported by González et al. (2012) with productions that oscillated between 200 and 6372 kg DM ha −1 year −1 in diverse rangelands with a wide range of variations in climate, livestock density and pasture improvements with fertilizations. Gómez Gutiérrez and Luis Calabuig (1992) studied several kinds of grasslands with scattered tree cover, determining annual productions lower than 500 kg DM ha −1 in many areas.
Plant growth depends on soil water availability that, in turn, is influenced by rainfall variations (Schnabel, 1997). Houérou and Hoste (1977) and González et al. (2012) found that the annual distribution as well as the interannual variations of precipitation had a significant influence in the correlation between precipitation and pasture production. The effect of rainfall variations on simulated pasture production for site 1 and site 2 are shown in Figs. 8 and 9, respectively. The graphs show annual pasture production over 300 years along with a 10 year window of results at the daily timescale that reflect the annual distribution of production. Annual pasture yield depended on annual rainfall amounts and the seasonal distribution, with periods of less yield corresponding to drier years, and greater productions in wetter years.  Figure 9. Simulated average pasture production and precipitation at Site 2, A) at the annual timescale for 300 years. B) for ten years at the daily timescale (the green shade represents +/-1 standard deviation of pasture production, and the blue bars is the rainfall) Fig. 9. Simulated average pasture production and precipitation at site 2; (A) at the annual timescale for 300 years; (B) for ten years at the daily timescale (the green shade represents ±1 standard deviation of pasture production, and the blue bars is the rainfall).
The seasonal distribution of rainfall did also influence pasture production. Accumulated antecedent precipitation before June was a good predictor of the yield regardless of the total annual precipitation. Years with low accumulated precipitation before June were less productive than years with higher accumulated precipitation (Table 6). For example, similar annual rainfall occurred in years 210 and 213; however, in the year 213 the rainfall of the last four months prior to June was higher, which resulted in a greater yield. In the year 215 a large amount of rainfall occurred after May, but pasture production that year was low.
Antecedent rainfall of the last 120 days before June was the variable that explained best the annual pasture production (r 2 = 0.73 and r 2 = 0.51, for site 1 and site 2, respectively). Shorter accumulation periods for antecedent precipitation had poorer correlations with yield, which can be explained because they are associated with less growing time and because as summer approaches there is an increase in evaporation losses.

Pasture production: spatial distribution
The spatial distribution of simulated pasture production varied greatly across the basin. Figure 10a presents the spatial  distribution of average production in the catchment over the entire 300 simulated years. Areas of higher production tended to have higher variability in their production ( Fig. 10b) as well as higher maximum and minimum productivities ( Fig. 10c and Fig. 10d). Productivity areas were persistent in time, with distributions determined by physiographic characteristics of the basin and the distribution of trees. A multiple regression analysis of pasture production with different variables showed that soil moisture, slopes, tree density, canopy cover, and upslope catchment area were the best predictors of production (r 2 = 0.81).
The distribution, composition and structure of plant communities are directly conditioned by spatiotemporal patterns in water availability (Asbjornsen et al., 2011) which is strongly determined by topography. In the study catchment the spatial distribution of the natural pastures was clearly influenced by the distribution of soil moisture. Areas with higher water availability had greater yield (Fig. 11a). Low yields were obtained if average soil moisture was lower than 0.150 m 3 m −3 . Slope also played a strong role in the distribution of yield. Topographically, valley bottoms and flat areas of the catchment were characterized by higher pasture production. Production decreased rapidly as slope increased Table 6. Annual pasture production at site 1 and site 2 (kg DM ha −1 ), annual rainfall (mm) and accumulated antecedent rainfall prior to 1 June (30, 60, 90, 120 days Figure 11. Scatterplot between average pasture production simulated and A) average soil moisture simulated, B) slope, and C) canopy cover. Fig. 11. Scatter plot between average pasture production simulated and (A) average soil moisture simulated, (B) slope, and (C) canopy cover.
( Fig. 11b). This is because in semiarid regions higher slopes are associated with reduced infiltration, enhanced drainage and production of overland flow (Cerdá et al., 1998). The importance of physiographic controls on soil moisture distribution and hence of pasture production in the study region was clearly documented in Ceballos-  and Van Schaik (2009), who demonstrated the importance of soils in low lying areas as water storages and the fundamentally different hydrologic regimes of hilltops, hillslopes, low areas and valley bottoms. Canopy cover exerted a strong control on pasture yield (Fig. 11c). An initial explanation is that pixels with high canopy coverage have higher interception of incident precipitation, more transpiration and therefore reduced soil moisture. This interpretation is, however, insufficient since the influence of trees on pasture production is a more complex issue that involves a number of processes not explicitly simulated in this study. For instance, trees may promote pasture production by enhancing soil fertility and structure or by providing a shaded and favorable microclimate. These factors were not explicitly simulated in this study. Still, it is known that in semiarid ecosystems, rainfall interception together with soil water uptake by trees in areas of high canopy cover would increase the competition for water resources between trees and pastures rather than enhance the productivity of pastures (Moreno, 2008). However, because the model used in this study does not incorporate many processes describing the overstory-pasture relationships such as the effect of vegetation on nutrients and on the soil microbial activity, we cannot conclude that tree canopy cover is strictly detrimental to the productivity of pastures. Indeed, several studies in the region show increased yield under trees as compared to open areas (Moreno, 2008). It has been observed that moderation of incident light could have a positive effect on crop productivity by altering the microclimate under trees, however this effect depends on antecedent conditions and the production of previous years (Gea-Izquierdo et al., 2009). Values of 13 % of canopy cover with 24 trees ha −1 were considered optimum for understory pasture production (Montero et al., 2008).

Climatic and physiographic factors
The degree to which the various controls discussed in the previous sections determine the distribution of pasture is not invariant. Precipitation is a main driver of total production ( Fig. 12a) in almost a linear fashion, but the spatial distribution of pasture is to a large extent controlled by topography, since the spatial variability of precipitation in the study area is very small. In Fig. 12a we distinguish between low, medium, and high production years. These years are clearly related to total precipitation amounts during the February-June period (50 to 150 mm of precipitation are associated with years of low production, 150 to 250 mm correspond to years of medium production and more than 250 mm yields high production). Rainfall is related to pasture growth through an associated increase in soil moisture available for uptake. While precipitation is related to production in a somewhat linear relationship, soil moisture is related to pasture productivity in a nonlinear, approximately sigmoidal relationship (Fig. 12b) that starts to reveal the effects of the heterogeneity of the terrain. Figure 12b suggests that the precipitation amounts only have a scaling effect on the relationship between soil moisture and pasture production. The functional form of this relationship or the ability of soil moisture to explain pasture production remains relatively unchanged.
Unlike rainfall, the distribution of soil moisture is affected by the heterogeneity of the terrain, but the strength of this effect is proportional to the amount of soil moisture, which is partially controlled by the amount of precipitation. For instance, low local slopes drive soil moisture by reducing flow velocity and by increasing the opportunity for infiltration; therefore, high production tends to be found in flatter areas of the terrain (Fig. 12c). The effect of the slope, though, is stronger during wetter years when soil moisture is higher and there is more opportunity for overland and subsurface redistribution of water. For drier years the ability of the local slope to explain the spatial variance of production decreases (Fig. 12c).
The relative position of a location in the drainage network, as defined by its upstream catchment area, is a non-local topographic control that also has a strong role in explaining the distribution of pasture production. More water is potentially drained at locations with a larger upstream catchment area, making them more prone to have a higher soil moisture content. Indeed, the productivity of a location increases with its upstream catchment area (Fig. 12d). Local drainage is defined by the small-scale topographic features of the surface that form a convergent network. During years of low precipitation, concentration of moisture in converging areas of the drainage network produces a very contrasting spatial distribution of pasture production. The strength of this topographic control during dry years can be assessed by its relatively high explanatory power of the total spatial variance of pasture production. For increasingly wetter years, the strength of this topographic control wanes and with it its explanatory power (Fig. 12d). The contribution of upstream inflows to total local soil moisture decreases as incident precipitation increases. This reduces the influence of the non-local topographic controls.
Overall, during years of abundant production of pasture the importance of upstream water inflows tend to be overwhelmed by relatively large inputs of precipitation. In these conditions local topographic controls such as low slopes that reduce local water drainage rates have a relatively higher influence in the observed pasture productivity. As precipitation inputs are reduced the importance of the lateral redistribution of water becomes more relevant and non-local controls such as the upstream drainage area becomes increasingly more explanatory of the distribution of pasture.

Conclusions
Ecohydrological spatially distributed models in conjunction with statistical weather generators are effective tools for simulating long-term pasture production dynamics and hydrologic conditions in semiarid rangelands, characterized by high spatial and temporal climatic and hydrologic variability. Results from this study contribute insight into the hydrologic and climatic controls that determine the spatial and temporal distribution of grasses and the expected range of pasture production in different areas at the watershed scale.
This study aims at informing rangeland management and promoting the sustainability of grasslands. Spatially, the general physiographic characteristics of the terrain are good predictors of pasture yield, but the distribution of the canopy overstory is also important. Valley bottoms and flat areas adjacent to slopes, which tend to have relatively high-soilmoisture contents, had the highest production in the study area. Tree canopy cover was found to be negatively related with pasture production, reflecting the importance of rainfall and light interception, as well as water consumption by trees, in the development of a grassy understory in semiarid rangelands.
The simulated pasture production in the study catchment ranged from 21 to 1030.9 kg ha −1 year −1 , which ranks it as a medium to low productivity compared to other Mediterranean rangelands. With the calculated yields, the introduction of supplemental fodder is necessary to maintain livestock. Although the interannual distribution of precipitation is a strong control on the variability of pasture yield, its seasonal distribution during the year is as important. Specifically, years with low rainfall from February to May showed limited yield even for years with relatively high annual precipitation.
The importance of topographic controls, as captured by the accumulated drainage area, becomes more relevant to explain the spatial distribution of pasture during years of low precipitation. This is because water inflows associated with lateral redistribution processes become a larger proportion of the total inflow into a location due to reduced precipitation inputs. The influence of lateral redistributions of water and therefore of the topographic structure of the watershed is reduced as spring precipitation inputs increase.
Although the model used in this study showed good performance in the simulation of water and vegetation dynamics in the study region and therefore provide confidence that the first order controls are captured, important processes, believed to play an important role in the long-term dynamics of pasture production, were not explicitly simulated. An example of these processes is the feedback between climatologic, ecohydrologic processes and the cycling of nutrients.