Plant’s gypsum affinity shapes responses to specific edaphic constraints without limiting responses to other general constraints

Harsh edaphic environments harbor species with different soil affinities. Plant’s responses to specific edaphic constraints may be compromised against responses to prevalent stresses shared with other semi-arid environments. We expect that species with high edaphic affinity may show traits to overcome harsh soil properties, while species with low affinity may respond to environmental constraints shared with arid environments. We quantified the edaphic affinity of 12 plant species co-occurring in gypsum outcrops and measured traits related to plant responses to specific gypsum constraints (rooting and water uptake depth, foliar accumulation of Ca, S and Mg), and traits related to common constraints of arid environments (water use efficiency, macronutrients foliar content). Plants in gypsum outcrops differed in their strategies to face edaphic limitations. A phylogenetic informed PCA segregated species based on their foliar Ca and S accumulation and greater water uptake depths, associated with plant responses to specific gypsum limitations. Species’ gypsum affinity explained this segregation, but traits related to water or nutrient use efficiency did not contribute substantially to this axis. Plant’s specializations to respond to specific edaphic constraints of gypsum soils do not limit their ability to deal with other non-specific environmental constraints.


Introduction
Harsh edaphic environments can be limiting for many organisms. As a result, the plant communities inhabiting these soils are characterized by sparse coverage and low biomass compared to those growing on more fertile soils in neighboring areas (Damschen et al. 2012;Escudero et al. 2015). Some plants living on stressful soils often have mechanisms to tolerate the toxicity imposed by certain elements (Moore et al. 2014), but other less stress-tolerant species can also colonize these environments without such specific mechanisms. This might result in differentiated strategies to deal with the harsh edaphic constraints for plant life found in these environments, potentially enhancing species coexistence and richness (Palacio et al. 2007;Escudero et al. 2015).
Plants adapted to harsh soils can be classified as edaphic endemics (hereafter specialists) or nonendemics (hereafter generalists). Specialists tend to show narrow edaphic tolerances, which restrict their ecological niche, while generalists have broader edaphic tolerances that allow them to survive in a wider array of soil types (Büchi and Vuilleumier 2014). It is commonly assumed that specialists have adapted to, and perform better, in environments with particularly stressful characteristics for plant growth than in other habitats (Levins 1968;Futuyma and Moreno 1988;Jasmin and Kassen 2007). However, some generalists can also thrive in these harsh habitats following an opportunistic strategy favored by environmental heterogeneity in space and time (Futuyma and Moreno 1988;Büchi and Vuilleumier 2014). Indeed, the coexistence of edaphic specialists and generalists is widely observed in harsh edaphic environments such as those derived from gypsum (Moore et al. 2014;Escudero et al. 2015), serpentine (Sianta and Kay 2019), granite (Murdy 1968) or dolomite (Mota et al. 2008).
Gypsum soils occupy over 100 million hectares worldwide (Verheye and Boyadgiev 1997). Gypsum ecosystems are mostly found in arid and semi-arid regions (Parsons 1976), limiting the establishment and survival of many plant species. Besides, gypsum also imposes other more specific edaphic stresses on plants, arising from its physicochemical properties. On the one hand, the low soil water and macronutrient (N, P, K) availability can be considered a common limitation that gypsum soils share with many other dryland environments. On the other hand, some of the particularly adverse physical limitations imposed by gypsum soils are the presence of a hard physical crust that limits plant establishment (Escudero et al. 2015) and its mechanical instability, high aggregation and low porosity (Bridges and Burnham 1980;Guerrero Campo et al. 1999a). These properties make gypsum a limiting substrate for vertical root penetration and development (Guerrero Campo et al. 1999b; Moore et al. 2014). Another adverse property of gypsum derives from its chemical composition (CaSO 4 ·2H 2 O), which generates an excess of Ca and S in the soil solution that can be detrimental for plant growth (Romão and Escudero 2005;Escudero et al. 2015). An excess of Ca in soil interferes with the uptake of other essential nutrients by plants due to Ca exchange with other soil ions (Guerrero Campo et al. 1999b), whereas S excess can be toxic for plants (Duvigneaud and Denaeyer-de Smet 1966;Ruiz et al. 2003).
In gypsum ecosystems, species with different degrees of gypsum affinity or specialization co-occur within the same plant community. These range from specialists only found on gypsum (gypsophytes) to a wide variety of generalists than can thrive on gypsum, but also on other lithologies (gypsovags). Plants living on gypsum exhibit different survival strategies that may respond to some of the harshest constraints of gypsum (e.g., high Ca and S concentrations or a hardphysical crust, high aggregation, presence of pure gypsum crystals and low porosity), or to other more general constraints shared with many arid ecosystems (e.g., low fertility and water availability). On the one hand, plant responses to deal with specific gypsum limitations could be related to facing chemical toxicity and soil physical resistance against root penetration and growth. An avoidance strategy to prevent chemical toxicity is the accumulation of Ca and S in plant tissues in response to their high concentrations (Ruiz et al. 2003;Palacio et al. 2014a). On the other hand, plants capable of overcoming rooting difficulties gain access to deeper soil layers with usually greater water storage during drought periods and lower inter-plant competition (Ryel et al. 2008). Plants living on gypsum can also show strategies to respond to other more common limitations, which could also be beneficial in other nutrient-poor and dry environments, such as an efficient nutrient acquisition or efficient water use.
Trade-offs among plant traits may emerge due to physiological constraints that limit the functional diversity of plant species. Trade-offs have been reported, for instance, between rooting depth, transpiration and water use efficiency (Brooks et al. 1997;Moreno-Gutiérrez et al. 2012). Plants living on gypsum may develop contrasting but equally successful strategies to cope with the stressful conditions imposed by the soils' physicochemical properties. Therefore, plants that safely accumulate excess ions (Ca and S), avoiding toxicity, might show a reduced ability to assimilate other essential nutrients such as N, P or K (Marschner 2012).
This study assesses whether there are trade-offs among traits so that plant's investment to face specific edaphic constraints is compromised against dealing with more prevalent stresses shared with other semi-arid environments. We hypothesize that a functional specialization to deal with specific gypsum constraints (e.g., deeper rooting and water uptake depth, Ca-S-Mg accumulation) may prevent water use efficiency and nutrient acquisition (e.g., higher transpiration and lower water use efficiency, lower N-P-K and C contents) due to the expected trade-off between the plant's investment in strategies to face specific and general constraints in semi-arid gypsum ecosystems.

Study area
We performed the study in a semi-arid Mediterranean ecosystem on gypsum soils located in the Vinalopó valley in southeastern Spain (Alicante, 38º 29 '39'' N; 0º 47' 00'' W). We selected flat areas to avoid topographical heterogeneity, demarcated within a radius of 13 km between 412 and 490 m.a.s.l. The dominant soil type was Keuper gypsum appearing abruptly in the form of intrusive outcrops, surrounded by other lithologies consisting mainly of limestone, but also clay and marl. Climate is semi-arid with an average temperature of 16ºC and a mean annual precipitation of 395 mm. Precipitation is strongly seasonal and falls mainly in spring (March-May) and autumn (September-November), with very low, or absent, precipitation in summer (June-August).

Evaluation of gypsum affinity and experimental design
We focused on 12 plant species commonly found on gypsum outcrops with a varying degree of gypsum affinity, including a wide phylogenetic diversity (Families in Table 1). For measuring gypsum affinity (i.e., gypsophily), we selected four localities in the same region where the boundary between the gypsum soil and the surrounding lithology (hereafter non-gypsum) was clearly demarcated. In each locality, we selected two contiguous subareas of approximately 1 ha, one within gypsum soil and another in non-gypsum soil (mainly limestone). Both types of substrates were closely located (< 100 m) in the four localities, sharing similar climatic conditions. We selected gypsum and non-gypsum areas to be as similar as possible in topography, avoiding areas with steep slopes. Sampling comprised 80 plots (150 × 150 cm) in each locality, except in one nongypsum locality with 79 plots. The plots were semirandomly distributed to occupy the 1 ha extension. The localities were sampled in four days periods twice per month between April 2019 and February 2020. Inside each plot, we identified all adult plants of the 12 target species (11,453 individuals) and measured each individual´s coverage employing the ellipse equation: Being a the semi-major diameter and b the semiminor diameter. Then, separately for each location, we calculated each species gypsum affinity (g) as the proportion of plant coverage found in gypsum as follow: Being Cg the coverage in gypsum areas and Cn the coverage in non-gypsum areas. The plants' coverage was estimated considering only the plants living alone, thus avoiding possible effect derived from the interactions between co-occurring plants not related to soil affinity. Gypsum affinity (g) values range from 0 to 1, where 0 indicates species found in the non-gypsum areas that never occur on gypsum and 1 indicates gypsophytes that only occur on gypsum. Species´gypsum affinity was determined as the mean g value for the target species in the four localities. This index gives a reliable measure of the degree of gypsum affinity for our studied community since it was estimated from in situ data. Finally, we measured traits in a total of 57 plant individuals of 12 species encompassing a wide gypsum affinity gradient (Table 1).

Plant water sources
Assessing rooting depth in the field can be challenging, but the analysis of xylem water's isotopic composition allows an indirect assessment of water uptake depth in woody plants (Dawson et al. 2002). In seasonally dry areas like the Mediterranean region, the intense evaporation of water from the soil surface during the hot dry summer produces isotopic fractionation resulting in an enrichment of the heavier oxygen ( 18 O) and hydrogen ( 2 H) isotopes in topsoil water. This evaporative isotopic enrichment generates a steep gradient in soil water isotopic composition with depth, with more enriched water in shallow soil layers and progressively less enriched water with depth (Allison et al. 1983;Smith et al. 1997;Teixeira et al. 2003). Little isotopic fractionation occurs during plant water uptake, especially for oxygen Ehleringer and Dawson 1992;Brunel et al. 1995;Dawson et al. 2002, but see: Ellsworth and Williams 2007; Barbeta et al. 2019) so the xylem water isotopic composition matches the mean isotopic composition of the different water sources taken up by active roots along the soil profile (Teixeira et al. 2003).
We analyzed the xylem water oxygen (δ 18 O) and deuterium isotopic composition (δ 2 H) of each plant in peak summer. We harvested lignified stem samples on August 14, 2017 early in the morning (7-9 am, solar time), once the plant is photosynthetically active but evaporative demand is low, to minimize stem water evaporation. The bark and phloem were scraped off the stems with a knife to avoid xylem water contamination with phloem water and organic compounds present in living cells and/or the bark (Ehleringer and Dawson 1992). After cutting, samples were immediately stored in individual airtight capped crystal vials and kept refrigerated in the field in a cooler until transportation to the lab where they were kept frozen at -80 ºC until extraction. Both xylem water extraction and stable isotope analysis of water were conducted at the Serveis Científico-Tècnics of the University of Lleida (Spain). Xylem water was extracted by cryogenic vacuum distillation (Ehleringer and Osmond 1989;Martín-Gómez et al. 2015). Sample vials were placed in a heated silicone oil bath (110-120°C), and connected with Ultra-Torr unions (Swagelok Co., Solon, OH, USA) to a vacuum system (ca.. 10 − 2 mbar) including U-shaped water traps in series that were cooled with liquid N 2 . The extraction time was 90 min. Captured water was then transferred into cap-crimp 2-ml vials, and stored at 4°C until analysis. The hydrogen and oxygen isotopic composition of the extracted xylem water samples was analyzed by isotope ratio infrared spectroscopy (IRIS) on a wavelength scanned cavity ring-down spectrometer (WS-CRDS) model L2120-i coupled to an A0211 Finally, we calculated the deuterium-excess (d-excess) for each xylem water sample using the relationship proposed by Dansgaard (1964).
Given that d-excess is derived from the relationship between δ 2 H and δ 18 O, it provides a precise measure to detect evaporative isotopic fractionation, and hence, differences in soil water uptake depth among plants.
Here, we assumed that low (more negative) values of d-excess imply enrichment in heavy isotopes, and thus plant utilization of intensely evaporated water from shallow soil layers (Allison et al. 1983).

Plant water use efficiency
We measured foliar δ 13 C and δ 18 O to infer the timeintegrated water use efficiency and stomatal conductance over the growing season in the studied plants. The carbon isotopic composition (δ 13 C) of the leaf is used as a timeintegrated proxy for intrinsic water use efficiency. The ratio between carbon uptake and stomatal conductance, i.e., the intrinsic water use efficiency (WUEi = A/g s ), can be estimated by the carbon isotopic fractionation occurring during CO 2 diffusion between the atmosphere and the sites of carboxylation, and during carboxylation itself (Farquhar and Richards 1984) The oxygen isotopic composition (δ 18 O) of foliar tissues provides a timeintegrated measure of stomatal conductance and, thus, cumulative transpiration (Barbour et al. 2000;Barbour 2007), being the foliar δ 18 O negatively correlated with transpiration (Farquhar et al. 2007). Foliar δ 18 O is unaffected by changes in photosynthetic rates (Scheidegger et al. 2000;Ramírez et al. 2009), but it is affected and included the water source isotopic signal (Sarris et al. 2013; Barbeta and Peñuelas 2017). When both carbon and oxygen isotopes are considered together, it is possible to separate the independent effects of carbon fixation and stomatal conductance on water use efficiency. Finally, it is important to remark that the transpiration rate is positively correlated with water uptake (Aston and Lawlor 1979;Cienciala et al. 1994).
In summer 2015, we collected 5 g of fully developed leaves from each plant individual, which were dried at 50ºC for 3 days and ground to a fine powder. We encapsulated 4 mg of ground leaf material into tin capsules for carbon isotope analysis (δ 13 C) and 0.2 mg into silver capsules for oxygen isotope analyses (δ 18 O). Samples were analyzed at the Centre for Stable Isotope Biogeochemistry, University of California, Berkeley (USA). Leaf δ 13 C was analyzed using an elemental analyzer (Carlo-Erba NS-1500, Milan, Italy) coupled to an isotope ratio mass spectrometer (Isoprime100, Elementar, UK). Leaf δ 13 C is expressed in delta notation (‰) relative to the Vienna Pee Dee Belemnite standard (V-PDB). Leaf δ 18 O was determined using an isotope ratio mass spectrometer (IRMS, ANCA/SL elemental analyzer) coupled to a Finnigan MAT Delta PlusXL IRMS Elemental Analyzer (Finnigan MAT, Bremen, Germany). Leaf δ 18 O is expressed in delta notation (‰) relative to the Vienna Standard Mean Ocean Water for δ 18 O. Long-term (3 + years) external precisions for δ 13 C and δ 18 O measurements of leaf material are 0.10 and 0.20 ‰, respectively.

Nutrient concentration in leaves
We measured the concentrations of all macronutrients, including those found in excess in gypsum (Ca, S and Mg) and those that can be limiting in gypsum and other semi-arid environments worldwide (N, P and K). We also measured the C concentration to assess differences in foliar stoichiometry due to the accumulation of certain ions. Leaves were dried at 50°C, milled, and P, K, Ca, Mg and S concentrations were measured using inductively coupled plasma optical emission spectrometry (ICP-OES, Thermo Elemental Iris Intrepid II XDL, Franklin, MA, USA) after microwave-assisted digestion with HNO 2 :H 2 O 2 (4:1, v:v). Foliar C and N concentrations were measured in an ANCA/SL elemental analyzer. Nutrient concentrations were measured at the Ionomic Service of CEBAS-CSIC (Murcia, Spain).

Phylogenetic relationships
All the statistical analyses considered the phylogenetic relationships among the studied plant species, as closely related species will tend to present similar traits and, therefore, should not be considered independent observations (Revell 2010). We assembled the phylogenetic relationships among the studied plant species with the R function "S.PhyloMaker" (Qian and Jin 2016), which matches a given species list (our plant community) with an expanded version of the time-calibrated angiosperm species-level mega-tree that includes more than 31,000 species with branch length representing chronological time (millions of years) (Zanne et al. 2014). Species not present in the mega-tree were randomly added to our phylogeny within their corresponding genera (scenario 3, described in Qian and Jin 2016). Finally, taxa not present in our community were pruned from our tree.

Statistical analyses
We used a multivariate approach to assess whether different plant strategies emerged using the measured variables. For this, we carried out a phylogenetically informed principal component analysis (herein, p-PCA), using all the measured variables (foliar Ca, Mg, S, N, P, K, C concentrations, d-excess of xylem water, δ 18 O leaf , and δ 13 C leaf ), including plant height as a variable in the p-PCA to account for possible effects derived from plant size. All variables were scaled previously to run the p-PCA with the "scale" R base function. The p-PCA was run using the R function "phyl.pca" in the R package "phytools 0.7.47" (Revell 2012). Finally, we conducted two phylogenetic generalized least square models (PGLS) using the first axis (PC1) and second axis (PC2) scores from the p-PCA as the response variable and gypsum affinity (g) as the predictor. PGLS is a comparative phylogenetic method that allows testing for the relationship between gypsum affinity and species strategy (defined by p-PCA axis), considering the expected covariance structure of residuals for a given phylogeny (our phylogenetic tree). The correlation structure was derived from a maximum likelihood estimate of Pagel's λ (Pagel 1997), using the "corPagel" function of the R package "ape 5.3" (Paradis et al. 2004). The PGLS was run using the "gls" function in the R package "nlme 3.1.147" (Pinheiro et al. 2019). All the analyses were performed using the statistical software R 4.0 (R Core Team 2019).

Results
Species differed widely in traits related to water uptake depth and foliar nutrients ( Table 2; phylogenetic relationships between the studied species are presented in Fig. 1). The first (PC1) and second (PC2) principal components of the p-PCA explained 43 % and 21 % of the total variance, respectively. Variables contributing the most to PC1 were foliar S, Mg, Ca concentrations and d-excess in xylem water (i.e., those specifically related to physical and chemical gypsum constraints), which showed highly negative loadings, and δ 18 O leaf , foliar C and, to a lesser extent, N concentration, which exhibited highly positive loadings ( Fig. 2; Table 3). Other variables such as foliar P and K concentration and δ 13 C leaf showed low absolute PC1 loadings ( Fig. 2; Table 3). The p-PCA also showed highly positive PC2 loadings for plant height, P and K concentration, and δ 1 3 C l e a f , and a negative PC2 loading for N concentration.
The PGLS analysis showed that the species scores along the PC1 of p-PCA were significantly and negatively correlated with gypsum affinity (standardized coefficient = -2.54 ± 0.64, F-value = 15.80, P-value = 0.003) (Fig. 3). Similar results were observed for individual relationships, with foliar Ca, S, Mg concentrations and d-excess of xylem water being positively correlated with g, whereas leaf δ 18 O and foliar C concentration were negatively correlated with g (see supplementary Table S1 for univariate responses). Results did not change substantially after excluding O. tridentata from the analysis (standardized coefficient = -1.90 ± 0.61, F-value = 10.00, Pvalue = 0.012), which indicates that the observed patterns were not exclusively driven by the extremely negative score of O. tridentata. Species with high gypsum affinity (low PC1 scores) exhibited strategies associated with traits having negative loadings, mainly related to high accumulation of Ca, Mg and S in leaves and acquisition of water from deeper soil layers. In contrast, species with low gypsum affinity measured in leaves are also presented (g 100 g − 1 ). * Species considered as gypsophytesDescription of studied shrub species, including gypsum affinity index (g), number of individuals of each species used to calculate g (Ng), number of individuals of each species used for traits measurement (Nt), and individual plant height (cm, mean ± SE) * Species considered as strict gypsophytes (high PC1 scores) showed strategies mainly defined by low cumulative transpiration (high δ 18 O leaf ), high foliar C and, to a lesser extent, high N concentration. This indicates that gypsum affinity (g values) . Dot colors range from dark purple (species with high gypsum affinity, g = 1) to light yellow (species with low gypsum affinity, g = 0). Arrows represent the loadings of each variable in the pPCA explained, at least in part, some of the variation along this PC1. On the contrary, we did not find a significant correlation between gypsum affinity and species scores along PC2 (standardized coefficient = -2.12 ± 1.18, F-value = 3.25, P-value 0.102), although foliar δ 13 C and, to a lesser extent, K concentration were negatively correlated to g when considering those variables individually (Table S1).

Main findings
Our results show that different strategies emerge to deal with the harsh edaphic environment imposed by gypsum. In this regard, the variation defined by the PC1 was mainly explained by the contrasting degrees of gypsum affinity of the target species. In one extreme of the PC1, the observed species strategy consists of responding to the edaphic constraints imposed by gypsum through deeper roots, hence overcoming the soil hardness, along with enhanced foliar Ca and S accumulation to deal with the soil chemical toxicity. The other extreme of this axis is defined by a combination of lower time-integrated transpiration and higher foliar C concentration and, to a lesser extent, a slightly higher N concentration. In agreement with our expectations, the lower scores of species with higher gypsum affinity on the PC1 indicate that their resource use strategy specifically responds to the edaphic constraints imposed by gypsum. However, contrary to our expectations, our results do not show traits related to plant responses to non-specific constraints (i.e., shared with other arid ecosystems) at the other extreme of the PC1 axis, such as high efficiency in water and nutrient use, although univariate analyses show that species with low gypsum affinity present high water use efficiency (foliar δ13C) and, to a lesser extent, K concentration (Table S1). Therefore, we conclude that species with high level of specialization respond specifically to the edaphic constraints imposed by gypsum, without hampering their response to prevalent constraints shared with other arid ecosystems.
Contrasting plant strategies depending on gypsum affinity Species with higher gypsum affinity may accumulate ions found in excess (S, Ca and Mg) as a mechanism to tolerate the high concentrations of these elements in gypsum soils or to adjust their osmotic potential to take up water from ionically extreme soils (Chen and Jiang 2010). This pattern is stronger for Ca and S but less consistent for Mg, as Mg accumulation ability is more species-dependent (Moore et al. 2014). Indeed, gypsophytes' ability to accumulate Ca, S and Mg ions has been previously demonstrated in Iberian gypsophytes (Duvigneaud and Denaeyer-de Smet 1966;Palacio et al. 2007;Cera et al. 2020), where this accumulation can occur in cell vacuoles directly in the form of gypsum crystals (CaSO 4 •2H 2 O) (Palacio et al. 2014). Our results suggest that the accumulation of inorganic S, Ca, and Mg may influence other physiological responses in plants living in this environment. On the one hand, the accumulation of inorganic elements can affect the foliar stoichiometry due to the high content of inorganic ions that may reduce, in turn, the foliar carbon concentration (Palacio et al. 2007). On the other hand, the accumulation of inorganic ions might help reduce the plant water potential, thereby improving soil water uptake (Flowers et al. 1977;Ajmal Khan et al. 2000). Moreover, deep soil layers usually remain wetter during long drought periods than shallow layers due to lower evapotranspiration. The greater access to water stored in deeper soil layers can also be associated with somewhat higher cumulative transpiration (lower δ 18 O leaf ) and Ca accumulation, although greater utilization of deep, non-enriched water may have also contributed to lower δ 18 O leaf values in species with higher gypsum affinity (Sarris et al. 2013). Contrary to our expectations, species without responses to gypsumspecific limitations do not show either a high nutrient or water use efficiency, despite being traits favorable to deal with common limitations in stressful dry environments. Instead, they seem to tolerate gypsum limitations without any specific strategies, showing a combination of low transpiration rate, potentially resulting from a low water availability derived from their limitations to access water in deep soil layers, and high foliar concentrations of C and, to a lesser extent N, potentially due to the reduced accumulation of excess elements such as S, Ca, and Mg.
Water source segregation based on gypsum affinity A far less explored topic is the potential vertical niche segregation regarding root scavenging for water at different depths in the soil profile, depending on the degree of species' gypsum affinity. It has been demonstrated that root systems typical of gypsovags face difficulties in penetrating gypsum soils (Bridges and Burnham 1980), while those of gypsophytes are better adapted to overcome gypsum structural difficulties, both at seedling (Romão and Escudero 2005) and adult stage (Palacio et al. 2014b). However, the traits or mechanisms that make specialists' roots better adapted to overcome gypsum physical constraints are still unknown. Our results suggest that species with different gypsum affinities have access to different water sources after considering their dimensions (height). Differential access to water pools can be considered a proxy for rooting depth by accounting for variation in species size (Schenk and Jackson 2002). These functional differences might segregate the water pool niches exploited by coexisting species depending on their gypsum affinity, thereby promoting the coexistence of individuals of species with different edaphic affinities on gypsum soils. Niche partitioning and complementary use of limiting resources reduces competition among coexisting plants and favors their coexistence (Chesson 2000), which may explain the final composition of the plant community on gypsum outcrops. Indeed, specialists and generalists coexistence is widely observed not only in gypsum ecosystems, but also in many other harsh edaphic environments such as serpentine (Sianta and Kay 2019), granite (Murdy 1968) or dolomite soils (Mota et al. 2008). Niche partitioning occurs in some of these systems, thereby stabilizing their high diversity, as observed in serpentines (Levine and HilleRisLambers 2009;Sianta and Kay 2019). However, the extent to which the coexistence of plants with contrasting degrees of edaphic affinity is due to niche partitioning must be further examined, not only in gypsum soils but also in other harsh edaphic environments.

Conclusions
Our study shows that individuals of species living on gypsum rely on different responses and strategies to deal with gypsum edaphic constraints based on their particular gypsum affinity. Species with high gypsum affinity rely on functional responses to deal with specific gypsum edaphic constraints (i.e., soil structural hardness and Ca and S excess). They respond to these edaphic limitations by accumulating Ca, S, and Mg highly abundant in gypsum soils and accessing water from deeper soil layers despite gypsum's strong physical constraints limiting root penetration and development. However, whether species with lower gypsum affinity rely on more generalist strategies such as higher water and nutrient use efficiencystrategies useful in other nongypsum arid ecosystemsremain uncertain.

Further research
Further research on edaphic generalists' physiological performance on gypsum soils will help understand the ecological filters that harsh edaphic environments impose on plants. However, our results do not show any compromise derived from edaphic specialization in terms of efficiency in water and nutrient acquisition and use. So the riddle of why specialists do not spread beyond their narrow edaphic optimum warrants further research by considering, for example, the importance of gypsum affinity on different fitness components, ranging from reproductive effort (traits related with flowering, fruit and seed production) to plant growth and survival. Reciprocal transplant experiments or greenhouse studies using gypsum and non-gypsum soils would be valuable for assessing specialists' performance in and off gypsum lithologies (Cera et al. 2020). It might also be interesting to explore whether the segregation of strategies observed between specialists and generalists to face the specific edaphic limitations imposed by gypsum can be generalized to other harsh edaphic environments, which may be fundamental to advance our understanding of plant species coexistence in these habitats.