Spatiotemporal variations in groundwater and evaporative demand drive ecophysiological functioning of a phreatophyte in drylands

Water is the main limiting factor for groundwater-dependent ecosystems (GDEs) in drylands. Predicted climate change (precipitation reductions and temperature increases) and anthropogenic activities such as groundwater drawdown jeopardize the structure and functioning of these ecosystems, presenting new challenges for their management. We developed a trait-based analysis to examine the spatiotemporal variability in the ecophysiology of Ziziphus lotus, a phreatophyte that dominates one of the few terrestrial GDEs of semiarid regions in Europe. We assessed morpho-functional and hydraulic traits along a naturally occurring gradient of depth-to-groundwater (DTGW, 2–25 m) in a coastal aquifer, and throughout the growing season of the species. Increasing DTGW and salinity negatively affected photosynthetic and transpiration rates, increasing plant water stress (lower predawn and midday water potential), and positively affected Huber value (sapwood cross-sectional area per leaf area), reducing leaf area and likely, plant hydraulic demand. However, higher atmospheric evaporative demand fostered higher transpiration rates and water stress. Differences in climatic conditions throughout the growing season drove temporal variability in Z. lotus responses along the DTGW gradient, with warmer and drier conditions promoting carbon assimilation and water loss more intensively at shallow water tables. This multiple-trait analysis allowed us to identify plant ecophysiological thresholds related to the increase in DTGW and evaporative demand during the growing season. These ndings highlight the existence of tipping points in the ecophysiological functioning of phreatophytic plants in drylands, which contribute to disentangle the functional responses of the related GDEs under groundwater detriment because of climate change effects.


Introduction
Water is an essential global resource for humans and ecosystems, particularly in arid regions where it is the most limiting factor (Newman et al., 2006). Arid and semiarid regions are characterized by low and shifting water availability across space and time (Eamus et al., 2013), thus vegetation has to live with water limitation or explore new water sources below ground (Arndt et al., 2001;Nardini et al., 2014). In this sense, groundwater reservoirs are crucial for the functioning of vegetation (O'Grady et al., 2006) in the ecosystems that have access to this hidden water source, the so-called groundwater-dependent ecosystems (GDEs) . GDEs of arid regions are highly vulnerable to alterations in the hydrological regime, because their structure and functioning depend on it . Groundwater condition, i.e., water quality and quantity, affects GDEs, and groundwater exploitation or pollution jeopardizes their structure and function as well as the species that constitute them (Zolfaghar et al., 2014;Eamus et al., 2015). How the structure and function of GDEs in arid and semiarid regions are affected by groundwater variations is a primary concern for scientists, managers, and policymakers who have to design sustainable plans to manage groundwater resources in the face of climate change (Klove et al., 2014).
Fluctuations in groundwater depth can be detrimental to the functioning of GDEs and the deep-rooted phreatophytic vegetation that tap groundwater (Naumburg, 2005). Groundwater drawdown can salinize both soils and water in arid regions due to the exclusion of salts by plants during water uptake or to the exposure of deeper and saltier groundwater (Jobbagy and Jackson, 2007; Runyan and D'Odorico, 2010).
Seawater intrusion, as an indirect effect of water table decline near the coast, is one of the main drivers of coastal aquifer salinization. Likewise, groundwater availability for plants can depend on salinity, which has shown substantial consequences in phreatophytic productivity, even inducing diebacks (Jolly et al., 1993;Doody and Overton, 2009; Runyan and D'Odorico, 2010). Even though salinity is a signi cant abiotic stress that intensi es drought impacts and water unavailability, there is little research on plant response to both groundwater salinity and depth (Kath et  Groundwater-dependent ecosystems are among the terrestrial ecosystems most vulnerable to climate change effects, and their ability to persist will depend on the resilience of phreatophytic vegetation to groundwater decline (Hultine et al., 2020). It is widely recognized that anthropogenic activities alter the groundwater regime, either directly through groundwater exploitation or indirectly through land-use change (Eamus et al., 2015(Eamus et al., , 2016, which in turn can promote soil and groundwater salinization (Jobbagy and Jackson, 2007; Nosetto et al., 2008). Additionally, future climate change, expressed in the Mediterranean basin by a reduction in precipitation and an increase in temperature (Giorgi and Lionello, 2008), will reduce groundwater recharge and raise evapotranspiration rates. Modelling carbon-water relationships will help us predict how hydrological changes can affect GDEs in terms of survival and productivity, thus addressing human impacts (Naumburg et al., 2005;Newman et al., 2006). To test vegetation response to altered water regimes, scientists usually resort to spatial gradients of aridity, altitude, water availability, and soil nutrients, among others (Lavorel and Garnier, 2002;Wright et al., 2004;Mitchell and O'Grady, 2015). Topography, for instance, can promote gradients in water availability, which cause critical variations in plant structure and function (Williams et al., 1996). The study of a species response to reduced water availability along environmental gradients will provide insight for identifying ecophysiological thresholds in phreatophytic vegetation Froend and Drake, 2006 (Cornelissen et al., 2003). In GDEs, plant functional traits are the vehicle to assess different aspects of ecosystem functioning as they respond to changes in the hydrologic regime . In this sense, an understanding of the connection between morpho-functional and hydraulic traits with groundwater characteristics (i.e., groundwater depth, conductivity, and temperature) will be crucial for predicting climate change effects upon GDEs.
Numerous morpho-functional traits such as Huber value (Hv), woody density, speci c leaf area (SLA), and gas-exchange rates show variation across depth-to-groundwater (DTGW) gradients in arid and semiarid  Zolfaghar et al., 2014;Osuna et al., 2015;Sommer et al., 2016;Nolan et al., 2017a). Hydraulic traits such as water potential are strongly correlated with DTGW gradients, as shown in phreatophytic oaks, eucalyptus, and acacias from California and Western and Central Australia (Carter and White, 2009;Osuna et al., 2015;Nolan et al., 2017a). Here, we explore a GDE dominated by the phreatophyte Ziziphus lotus (L.) Lam. (Rhamnaceae) in a small coastal plain in the southeast of Spain where spatiotemporal variations in groundwater salinity and temperature were also assessed. We evaluated the relationships among a broad suite of morpho-functional and hydraulic traits including stem water potential, gasexchange rate, intrinsic water-use e ciency (WUEi), Huber value (Hv), wood density, and speci c leaf area (SLA), across a naturally occurring DTGW gradient related to distance from the coastline. We also assumed that seawater intrusion could more adversely affect plants near the coast. Thus, we hypothesized that spatiotemporal uctuations of both groundwater availability and quality would drive differences in the ecophysiological functioning of Z. lotus. These differences could help us to identify ecophysiological thresholds, which will provide valuable insight to face upcoming management challenges in GDEs. To test these hypotheses, we address the following speci c questions: Are there spatiotemporal variations in plant morpho-functional and hydraulic traits? Do these variations respond to groundwater conditions?
Is there any discernible threshold in the ecophysiological functioning of Z. lotus? What factors drive the threshold? 2. Methodology

Site description
The study was conducted on a coastal plain at the western part of the Cabo de Gata-Níjar Natural Park, southeastern Spain (Fig. 1). The plain is dominated by the phreatophyte Z. lotus that thrives with shallowrooted Mediterranean shrubs such as Lycium intricatum, Salsola oppositifolia, and Withania frutescens (Tirado, 2009). Vegetation shows a dispersed pattern with patchy distribution typical of arid and semiarid Mediterranean regions, where Z. lotus is associated with biodiversity islands (Tirado, 2009). Z. lotus is responsible for most of the photosynthetic activity during summer, whereas the rest of the vegetation constituting the island grows in winter, entailing a replacement in the drivers of the primary productivity of the ecosystem (Guirado et al., 2018). The climate is characterized as Mediterranean and semiarid, with hot and dry summers and mild winters, and mean annual precipitation of 200 mm which is unevenly distributed during Spring and Autumn (Machado et al., 2011). The shallow aquifer upon which Z. lotus depends comprises Plio-Pleistocene conglomerate, with sandstone beneath it and Pliocene marine marl at the base. The geology results from the sedimentary basin of Sierra Alhamilla mountains (1000 m.a.s.l) and from marine deposits from the Quaternary period (Vallejos et al., 2018). Eight bores located along the study area form a net for groundwater observation that discerns between 3 sites (east plain, west plain, and the seasonal stream that crosses it) and shows a natural occurring DTGW gradient based on coastline distance and topography. Vegetation sampling was made on a total of 16 Z. lotus individuals selected next to each bore (two per bore) ( Fig. 1) in three speci c periods of 2019 growing season: late-Spring (May), mid-Summer (July), and late-Summer (September).

Hydrologic and climatic measurements
Each bore contained two sensors (Hobo U20 Water level logger and Hobo U24 conductivity logger, Onset Comp. Coorp., Bourne, MA, USA) to obtain DTGW, electrical conductivity (i.e., salinity), and groundwater temperature (T GW ) every 15 minutes since May 2019. For regression analysis, we obtained mean values from each of the sampling periods. In the same way, we collected climatic data from Almería airport meteorological station every day (AEMET station), although just monthly precipitation (P) and mean monthly temperature (T air ) were used.

Hydraulic and morpho-functional traits
We measured water potential during the growing season at predawn (Ψ pd ) and midday (Ψ md ) in four stems on each of the 16 individuals using a Scholander pressure chamber (SKPM1405, Skye Instruments, Powys, UK). Measurements were taken before sunrise for Ψ pd (from 06:00 to 07:00 hours in May and July and from 06:30 to 07:30 hours in September) and during the peak insolation for Ψ md (between 13:00 and 14:00 hours). Mean values for each plant and period were calculated, and the maximum daily range (ΔΨ max ) was derived afterward as the difference between Ψ pd and Ψ md . We measured leaf gas exchange in 8 leaves per plants, covering all variation within the canopy, between 10:00 and 13:00 hours on the same days as water potential was measured. A portable infrared gas analyser (Li-6400XT; LI-COR Inc., Lincoln, NE, USA) was used with the following conditions in the chamber to standardise all measures: ow rate, 400 µmol s − 1 ; CO 2 concentration, 400 µmol mol − 1 ; and light intensity, 1800 µmol m − 2 s − 1 .
Ambient temperature was kept, which varied between 25-30 ºC. We obtained photosynthetic rate (A), stomatal conductance (g s ), transpiration rate (E), vapour pressure de cit (VPD), and WUEi was calculated from the ratio between A and g s .
Finally, we cut three branches per plant in July from which all leaves were removed. We measured sapwood cross-sectional area with a digital calliper in the base of each branch and estimated wood density as the volume of a piece of branch divided by its dry weight (after 48 h at 60ºC). We scanned all the leaves with a digital leaf-area meter (WinDIAS, Cambridge, UK) to calculate total leaf area (LA) per branch and used ten of the leaves to estimate the SLA of the plants, which represents the relationship between the LA and its dry weight (after 48 h at 60ºC). We calculated the Hv per plant from the ratio between the mean sapwood cross-sectional area to the mean LA.

Data analysis
We applied a two-way ANOVA for each groundwater characteristic and functional trait to assess intraspeci c variability, both temporal (between sampling periods) and spatial (between sampling sites).
Because SLA, Hv, and wood density were only measured once, we performed a one-way ANOVA for these traits. All traits were log-transformed except for water potentials due to the negative nature of their values.
We undertook Tukey's HSD post-hoc test after signi cant differences were found. We performed multiple bivariate linear regressions to test whether a single regression could describe individual functioning.
Some regressions were made with mean values, as variability over time was not observed, whereas others were made with monthly data to detect seasonal patterns. Finally, we analysed multiple trait relationships across all variables with a principal component analysis (PCA). Traits were scaled prior to the analysis to obtain a unit variance. Spearman correlation analysis was applied, and the contribution of each trait in the PCA was assessed to select those variables that provide the best representation and improve the analysis. Because of that, SLA, WUEi, and wood density were not included in the nal analysis. We performed all analyses in R 3.5.2 (R Core Team 2018).

Spatiotemporal variations in groundwater
We observed signi cant differences in DTGW, salinity, and T GW between sites, across the growing season, and for their interaction (P < 0.001; df = 7, 4). These variables increased during the growing season, although with different patterns. First, DTGW that ranged from 2.1 m (bore 1) to 25.4 m (bore 8) (Fig. 2a) increased across the growing season, although not substantially (Online resource 1). It was just at the inner-plain sites where an average increase of 18 cm was observed at the end of the season (bore 8).
Near the coast, we observed more noticeable temporal uctuations although these did not entail overall DTGW increments (Fig. 2b, c, d, and Online resource 1). Second, T GW gradually increased during summer despite its narrow range (Online resource 1 and Online resource 2). These rises mainly affected bores with the shallowest water tables. bore 1 showed wider uctuations, whereas bore 2, with the lowest T GW , had the steepest increase. Finally, groundwater salinity, which ranged from 3360 µS/cm (bore 4) to 11000 µS/cm (bore 7), increased in bores 1, 3 and 7, but particularly in bore 7 where a rise in almost 1000 µS/cm was observed (Online resource 1 and Online resource 2). For these three groundwater properties (depth, temperature, and salinity), uctuations were larger near the coast where water tables were shallower (bores 1 to 3) than in the other bores (Online resource 2).

Spatiotemporal variations in plant morpho-functional and hydraulic traits and their relationship with groundwater
Morpho-functional and hydraulic traits also showed signi cant differences between sampling periods, sites, and the interaction between them (Table 1 and Online resource 3). Overall, gas exchange (A and E) in Z. lotus leaves was higher in summer (July and September) and at those sites with the shallowest water tables. Regarding water loss, bores 1 to 4 (DTWG < 11.6 m) showed the highest g s , especially during July and September when it reached 0.42 ± 0.03 mol H 2 O m − 2 s − 1 , whereas bores 5 to 8 (DTGW > 14.0 m), the lowest values. It is also noticeable that high rates of E for plants from bores 1, 2, and 3 occurred in July and September, but also from bore 8 (25.3 m). Nevertheless, values of A showed signi cant differences in summer just at some locations (interaction term, P < 0.001). Individuals next to bores 2 and 5 (with a DTGW of 7.3 and 14.0 m respectively) had higher photosynthetic rates in July, whereas plants near bores 6 and 8 (with 19.3 and 25.3 m respectively) showed lower values in summer (Table 1). In general, individuals next to bores 1 and 2 had the highest values of A (15.51 ± 1.12 and 24.88 ± 2.36 µmol CO 2 m − 2 s − 1 respectively), whereas bore 8 showed the lowest rates (7.16 ± 1.33 µmol CO 2 m − 2 s − 1 ). However, A did not show differences between months. Contrary to A, WUEi was low at not only the shallowest water tables, but also at the deepest ones. Regarding water potential, more negative values of both Ψ pd and Ψ md were observed in July and September at sites with the highest DTGW (bores 5 to 8). Ψ pd ranged between − 0.32 ± 0.02 MPa in May to -1.55 ± 0.09 MPa in September (at bore 2 and bore 8, respectively), whereas Ψ md showed values between − 1.18 ± 0.04 MPa in May to -3.13 ± 0.10 MPa in July (bore 4 and bore 8, respectively). Hv showed signi cant differences across sites (P = 0.027) (Online resource 4). Plants at bore 1 had the lowest value of Hv (3.58 ± 0.08), whereas plants at bores 7 and 8 had the highest ones (11.40 ± 0.22 and 9.34 ± 0.84 respectively). Neither SLA nor wood density showed signi cant spatial variability. Table 1 Summary of mean values of traits (± standard error) from plants next to each bore in the three sampling periods: May, July, and September. Depth-to-groundwater (DTGW) of each site is showed as well as the significant differences (P < 0.05) between months in each site (different letters). Photosynthetic rate (A), stomatal conductance (g s ), transpiration rate (E), intrinsic water-use efficiency (WUEi), predawn (Ψ pd ) and midday (Ψ md ) water potential, and vapour pressure deficit (VPD). Morpho-functional and hydraulic traits signi cantly responded to spatiotemporal variations. First, bivariate linear regressions revealed a weak negative relationship with DTGW for most gas-exchange traits during the growing season (Fig. 3), except for WUEi. By contrast, no relationship was observed between these traits and groundwater salinity (Online resource 5). Regarding hydraulic traits, Ψ pd was the only variable that showed a signi cant linear relationship to both DTGW and salinity, in which Ψ pd but not Ψ md was signi cantly lower when DTGW and salinity were large ( Fig. 4 and Online resource 6). Our results also showed that at large DTGW, plants had higher Hv values than when DTGW was small (Fig. 5), even though wood density and SLA did not respond to groundwater gradients (Online resource 7). Therefore, DTGW was the main variable related to spatial variation in most single morpho-functional and hydraulic traits.
Vapour pressure de cit, representing temporal variations in climatic conditions, showed a positive correlation with E during the growing season, particularly in the warm and dry summer. In May, both E and VPD showed lower values, with little variability across bores, whereas in summer (July and September), the increase in VPD was parallel to the rise in E (Fig. 6a). The general increase in VPD during the season enhanced transpiration rates more over the shallowest water tables than at the deepest ones (Fig. 6b). However, VPD did not show any signi cant relationship with other traits related to gas exchange (Online resource 8). The overall increment of VPD from spring to summer was related to more negative Ψ pd and Ψ md values, as shown in the regression analysis ( Fig. 7a and b).
Temporal analysis of the relationships between traits also revealed that A, E, and g s were positively related to each other, particularly during spring, whereas WUEi (= A / g s ) was positively related to A and negatively related to g s in summer exclusively (Online resource 9). Our results also showed a negative relationship of Ψ pd with these gas-exchange traits both in spring and summer. As water availability decreased (lower Ψ pd ), A, g s , and E were reduced, but no response was observed with an increase of plant stress (lower Ψ md ) at any time (Online resource 9).

Multiple trait relationship for identifying ecophysiological thresholds
PCA revealed multiple trait relationships that were not identi ed with simple regression analysis. The two rst components of the PCA explained 63.5% of the variation across plants (Fig. 8a). The rst component (PC1), accounting for 37.6% of the total variation, showed strong loadings for climatic variables (i.e., T air , precipitation) as well as hydraulic traits (i.e., Ψ pd and Ψ md ) and E. The second component (PC2) explained 25.9 % of the variance and showed strong loadings for groundwater traits (particularly DTGW but also groundwater conductivity), A and g s . Groundwater salinity and Hv also contributed to PC2, although to a lesser extent. As a result, axis 1 showed a temporal gradient from the warmest and driest months that overlap to each other (July and September) with higher E and VPD, to the mild and humid spring (May), when water availability was higher (high Ψ pd ) and plant stress lower (high Ψ md ) (Fig. 8b). By contrast, axis 2 showed a DTGW gradient ( Fig. 8c and d) where plants closer to the water table exhibited higher A and g s but lower Hv. The PCA revealed two distinct clusters based on groundwater characteristics (DTGW, salinity) and their associated gas-exchange traits (A, g s ): one for plants at sites with shallow DTGW (< 12 m, Fig. 8c), and the other for plants at sites with deep DTGW (> 14 m, Fig. 8d).

Discussion
In this study, we examined the ecophysiological response of the long-lived phreatophyte Ziziphus lotus to a DTGW gradient, in a coastal GDE of the Mediterranean basin. We found that DTGW and salinity had a signi cant effect on the ecophysiological function of this phreatophyte, as hypothesised. We further found that some traits were more strongly correlated to uctuations in DTGW and salinity (e.g., A and g s ), whereas others were more strongly related to seasonal uctuations in climatic conditions (e.g., E, Ψ pd , Ψ md ). By applying a multiple-trait approach, we were able to identify plant ecophysiological thresholds related to the groundwater characteristics and seasonality throughout the growing season.

Spatiotemporal variations in Z. lotus' morpho-functional and hydraulic traits and their relationship with groundwater
Our ndings revealed spatiotemporal variations in Z. lotus morpho-functional and hydraulic traits, which were related to both groundwater and seasonal climatic conditions. The spatial variability in DTGW might explain the response patterns of gas exchange throughout the growing season. Increasing DTGW negatively affected carbon assimilation and water loss, as has been previously observed in a variety of GDEs ( lotus transpiration rate did not decline in summer; in fact, it increased with VPD, more signi cantly at shallow water tables, suggesting that summer conditions could induce higher rates when su cient groundwater is available Eamus and Prior, 2001), and that groundwater availability to the plant depends on climatic conditions. Despite the risk of hydraulic failure due to this anisohydric behaviour (Torres-García et al., 2021) and the physiological limitations of tapping water from deep sources, Z. lotus plants can maintain high gas exchange under current conditions.
The naturally occurring DTGW gradient also explained the spatial variability in Ψ pd and responses to differences in water availability. Ψ pd largely re ects the water potential of the rooting area (Hinckley et al., 1978) and indicates groundwater access by plants when values are relatively high (Carter and White, 2009). Although Z. lotus plants showed values higher than − 1.55 MPa, which is high given the solute potential, we found a negative trend of Ψ pd not only with increasing DTGW but also in salinity.
Groundwater electrical conductivity increased with DTGW away from the coast, which could be due to a marine incursion during the Holocene that penetrated the inner parts of the plain, constituting a lagoon which dried up over time and increased the salinity of the area (Vallejos et al., 2018). This result is contrary to our assumption that seawater intrusion could affect salinity near the coast. Instead, we found that the combination of deep groundwater and high salinity away from the coast might promote water stress in the root zone (Kath et al., 2015). Physiological and even anatomical constraints to accessing deep and salty groundwater might have induced lower Ψ pd and reduced water uptake despite Z. lotus continuous groundwater use even at these locations. Therefore, it is not a recent process of seawater intrusion that induced differences in Z. lotus population, but a past event that fostered different salinity conditions across the landscape. where daily and seasonal groundwater uctuations are minor, phreatophytes run the risk of maximizing productivity over safety (Hultine et al., 2020), which can also be fostered by the observed anisohydric behaviour. Continued transpiration under VPD higher than 4 kPa might lead to hydraulic failure if water becomes scarce (Drake et al., 2018), which is more likely with deeper DTGW. Therefore, being an anisohydric phreatophyte in arid and semiarid regions seems to be a risky option, which can only be overcome by the resilience capacity of the plants (Hultine et al., 2020).
Different responses observed in Z. lotus transpiration rates could also be generated by differences in xylem traits (Attia et al., 2015). Our results revealed that Hv (the ratio of sapwood cross-sectional area to leaf area) was higher at deeper groundwater sites, as already reported for other phreatophytes of mesic (Zolfaghar et al., 2014) and xeric environments (Carter and White, 2009). Thus, variations in Hv with the DTGW gradient would derive from differences in leaf area instead of wood density, suggesting that plants with less reliable groundwater supply (deep DTGW) make smaller investments in leaf area than plants at shallow sites (Zolfaghar et al., 2014). This mechanism might allow Z. lotus to cope with reduced water availability by decreasing their hydraulic demand and, therefore, their transpiration rates at a canopy level (Gazal et al., 2006;Carter and White, 2009;Zolfaghar et al., 2014). Indeed, reductions of aboveground biomass are acknowledged to be a common adaptation when plants cannot overcome the anatomical and functional adaptation cost of water scarcity (Naumburg, 2005). In contrast to Hv, wood density was largely independent of groundwater because it depends on development of modi ed cell types (e.g., xylem vessels, bres) (Lachenbruch and McCulloh, 2014). Likewise, our results showed that SLA was independent of DTGW, as has been re ected in some studies along water availability gradients (Nolan et al., 2017a). In this case, as SLA refers to the ratio of leaf area to leaf dry mass, or the inverse of leaf thickness (

Ecophysiological thresholds and future considerations
Assessing the expression of multiple traits provides tools to predict patterns of change in GDEs in response to variability in groundwater and across seasons (Hultine et al., 2020). A multiple-trait analysis revealed that the variability observed in the functioning of Z. lotus could be explained by the combination of both temporal variations in climatic conditions during the growing season of the species and the spatial differences in groundwater characteristics of the study area. Temporal differences from spring to summer showed a decrease in water potential with increased transpiration rates, promoted by environmental conditions (lower humidity, higher temperatures, and evaporative demand). This response could have fostered evaporative cooling, regulating leaf temperature for maintaining the plant carbon balance (Drake et al., 2018) and suggesting the decline in water potential was insu cient to indicate water stress. Thus, Z. lotus plants could avoid extreme thermal stress that can damage the photosynthetic machinery whilst preventing a steep decline in photosynthetic rate. However, su cient water availability is required to maintain evaporative cooling, which is essential under ongoing increases of both mean air temperatures and the severity of heat waves (Urban et al., 2017).
By contrast to temporal uctuations in ecophysiology, the ecophysiological functioning of Z. lotus across space was explained by the combination of groundwater availability (mainly determined by DTGW) and salinity (expressed by electrical conductivity). Salinity is commonly present in arid ecosystems with phreatophytic vegetation because of reduced precipitation, which prevents leaching of salts, and evaporation, which leaves salts behind (Glen et al., 2013). We found that the DTGW gradient coincided with a salinity gradient such that the deepest groundwater was also saltiest. Without the ability to discriminate between these characteristics, we observed that higher groundwater salinity combined with larger DTGW affected the ecophysiology of Z. lotus and promoted remarkable differences along the naturally occurring gradient. We identi ed a response threshold at 12-14 m, mainly promoted by differences in gas-exchange rates, which is consistent with previous studies about the species (Torres-García et  . The concern is also whether a depletion in groundwater level would exceed the root growth rate (Orellana et al., 2012), or even if temporal uctuations would have a long-term impact on plant ecophysiology. In this sense, phreatophytes that obtain groundwater from deep water tables and that already experience physiological constraints (e.g., over 14 m in the case of Z. lotus), could be intensively jeopardized by groundwater variations in the future.

Conclusions
In this research, we assessed spatiotemporal variations both in groundwater properties of a GDE in a semiarid region and in the morpho-functional and hydraulic traits of the phreatophyte that dominates this ecosystem: Ziziphus lotus. The naturally occurring DTGW gradient and associated monitoring eld station have provided an interesting scenario to assess ecophysiological differences related to water availability for phreatophytic vegetation. Here, we show that both groundwater depth and salinity are highly connected to the ecophysiological functioning of phreatophytic vegetation in drylands.
Nevertheless, no evidence of seawater intrusion seemed to affect Z. lotus plants, and groundwater salinity could be related to past events of seawater rise. Differences in climatic conditions throughout the growing season drove temporal variability in Z. lotus response, with summer conditions promoting carbon assimilation and water loss in this winter deciduous phreatophyte, more intensively at shallow water tables. The multiple-trait analysis led to identifying spatial and temporal ecophysiological thresholds that depend on both groundwater availability and atmospheric evaporative demand. Under the expected reductions in groundwater reservoirs as consequence of both climate aridi cation, and the increase in groundwater consumption and drawdown by human overexploitation, understanding the structure and functioning of GDEs of arid and semiarid regions and de ning ecophysiological thresholds of their phreatophytic vegetation will provide valuable insight to face upcoming management challenges.
Declarations Figure 1 Distribution of the boreholes (squares 1 to 8) and the related plants of Ziziphus lotus (circles, n = 16) on the coastal plain of Cabo de Gata-Níjar Natural Park, southeastern Spain.    and R2, the goodness of the t. Signi cance of the regression: **P < 0.01.  Bivariate linear regression between vapour pressure de cit (VPD) and predawn (a) and midday (b) water potential (Ψpd and Ψpd, respectively). Monthly values per plant are displayed ± standard error. Colours and shapes represent sampling periods: May (Green triangles), July (Yellow circles), and September (Red squares). R2 represents the goodness of the t. Signi cance of the regression: ***P < 0.001.