The Mediterranean forests studied are located in Sierra Morena Mountains, in the south of Spain (Córdoba province). The area is characterized by a continental-Mediterranean climate with cold, wet winters and dry, warm summers. Mean annual temperature is 17.6°C (with maximum values in summer reaching 40°C) and mean annual precipitation is 536 mm (with a 3-month period in summer without rainfall; data from AEMET for the years 1971–2000¹). Several shrub and arborescent species, such as Cistus albidus and Quercus coccifera, are abundant in drier soils, while broad-leaf deciduous trees, such as Alnus glutinosa and Fraxinus angustifolia, are dominant in moister soils (see Supplementary Table 1 for details of the 45 studied species). One cultivate species (Cydonia oblonga) has been included, because is naturalized in the study area and could potentially alter the biogeochemical niche of the coexistent species. Twelve sampling sites distributed over four different south-facing slopes were selected along a topographic gradient (from ridges to valley bottoms; see Supplementary Figure 1) with the aim of spanning a broad range of variations in soil resource availability, mainly of soil water (Supplementary Figures 2, 3 and de la Riva et al., 2016b). The species composition was recorded by measuring the cover of each woody species intercepted by four 20-m transects in each sampling site.

For trait measurements, we selected all the species appearing in the sampling transects, excluding only those with a relative abundance below 1% (in these cases it was difficult to find at least six individuals per species in the sampling site). This gave a total of 45 selected species, many of them appearing in more than one sampling site (Supplementary Table 1). We measured different leaf traits associated with resource uptake and conservation: two key functional traits related to morphology and associated with light capture and growth rate [specific leaf area (SLA)] and stress tolerance [leaf dry matter content (LDMC)], all the macronutrients (N, C, P, K, S, Ca, and Mg) and the isotopic ratio of carbon, which is related with gas exchange and water-use efficiency (Pérez-Harguindeguy et al., 2013).

Six samples of leaves from different individuals per species and site were collected in spring. Specific leaf area (leaf area per unit dry leaf mass; m² kg^-1) and LDMC (dry mass per unit of water-saturated fresh mass; g g^-1) were measured according to the methods recommended by Pérez-Harguindeguy et al. (2013). Chemical composition was determined for a mixture of leaves, from six different individuals per species and sampling site (except N and C which were measured in each individual). N and C concentrations were measured using an elemental analyser (Eurovector EA, 3000; EuroVector SpA, Milan, Italy). The macronutrients P, K, S, Ca, and Mg were extracted by wet oxidation with concentrated HNO₃ under pressure in a microwave digester, and analyzed by ICP-OES (Domínguez et al., 2012). Carbon isotope ratio (δ¹³C; ‰) was measured by combustion at 1020°C using a continuous flow isotope-ratio mass spectrometry system by means of Flash HT Plus elemental analyzer coupled to a Delta-V Advantage isotope ratio mass spectrometer via a CONFLO IV interface (Thermo Fisher Scientific, Bremen, Germany). The analytical measurement errors were ±0.1‰ for δ¹³C.

For all the statistical analyses, plant species were sorted into different groups according to (I) their topographical position (environment): Ridge Forest (hereafter RF), Middle-slope Forest (hereafter MF) and Riparian Forest (hereafter RiF); (II) leaf habit: Evergreen (hereafter Ev), Deciduous (hereafter De) and Semideciduous (hereafter Sd); and (III) growth form: Shrubs (hereafter Sr), Trees (hereafter T), Arborescent-Shrubs (hereafter ST), and Climber (C) sensu lato. Semideciduous species are considered those that have the capacity to drop their leaves under severe drought conditions (Zunzunegui et al., 2005; Ciccarelli et al., 2016), and therefore they have functional differences with respect evergreens.

A general principal components analysis (PCA) was performed with the whole set of leaf structural and chemical traits (10 leaf variables) for the 100 observations (measurements) made of woody plants (of 45 species, as some species appear in more than one sampling site), to study the degree of leaf trait variation among them. We used linear mixed models (considering species as a random variable) to assess whether the PCA scores of the first, second and third components (see details in Supplementary Table 2) differed among all possible combinations of the following fixed explanatory factors: environment, leaf habit, and growth form.

The total niche space of the community was calculated by the estimation of the n-dimensional hypervolume (Blonder et al., 2014), from the trait space occupied by the total species group in the different categories (environment, leaf habit, and growth form; the number of observations for climbers was not enough to calculate the hypervolume). In order to reduce the number of dimensions (which is recommended for this analysis), we used the first three PCA axes to calculate the hypervolume for each group, using a multidimensional kernel density estimation (KDE) procedure (see Blonder et al., 2014 for mathematical details). The units of the hypervolumes are reported as the standard deviations of centered and scaled log-transformed trait values, raised to the power of the number of trait dimensions used (sd^number of dimensions). We also calculated the overlap between the hypervolumes of each group with the correlation analysis of the “hypervolume” package, which compares the similarity between different hypervolumes using the Sørensen index (see Blonder et al., 2015). A rarefaction analysis was performed to control for the effects of species richness on the hypervolume. Thus, for each environment type, leaf habit and growth form, we built 100 randomized communities composed of species drawn (nine observations) from the species pool of that group. Then, we calculated the hypervolume of each sample and performed a one-way ANOVA to compare the hypervolumes of the groups, independently of species richness. In addition, to assess the functional trait overlap within and among groups (leaf habit and growth form) within each environment type, we calculated the hypervolumes and the Sørensen index between each pair of the most representative groups along the gradient.

All these analyses were conducted in the R 2⋅10⋅0 statistical platform (R Development Core Team, 2011), using the packages “vegan” (Oksanen et al., 2013), “nlme” (Pinheiro et al., 2015), and “hypervolume” (Blonder et al., 2014).

The overall PCA showed a clear separation among the 45 woody species in the volume defined by PCA₁ (accounting for 48.2% variance), PCA₂ (12.6% variance), and PCA₃ (8.5% variance) (Figure 1). The PCA₁ mostly reflected a gradient of increasing LDMC and C, and decreasing SLA, K and Mg; the PCA₂ was mostly related in one extreme (negative values) with high N, P, and S and at the opposite extreme (positive values) with high values of Ca and Mg; the PCA₃ was represented mainly by variations in N, S, P, and K (see variables scores in Supplementary Table 2).

The woody species inhabiting the ridge forest (RF) had different trait values (well-separated according the first and the third PCA axes) from those living at the Riparian forest (RiF) (Figure 2). There was also a separation in the trait PCA between species grouped by leaf habit: Deciduous (De) were different from Sd and Ev. With regards to growth form, there were differences among the arborescent shrubs (for PCA₁), and trees (for PCA₃) with arborescent shrubs and climbers (Figure 2); while no significant differences were found for PCA₂.

The results from the n-dimensional hypervolume approach showed that the functional space was greatest for the riparian communities from riparian forest (Figure 3A), the group of deciduous species (Figure 3B), and the growth form “trees” (Figure 3C), in the three plant dimensions (PCA axes 1, 2, and 3). In addition, after standardizing for species richness, the functional space showed significant variation along the topographic gradient: the hypervolume was significantly greater at the riparian forest (P < 0.001; Figure 4). Semideciduous species showed the lowest hypervolume, while deciduous species had the highest. Among the growth form types, trees had the highest hypervolume (Figure 4).

The degree of overlapping among the hypervolumes of the different communities was variable, ranging from 28 to 90% (Figure 3D). Attending to the environment, the overlap degree of the ridge forest RF was greater with middle MF than with Riparian forest RiF. Semideciduous type showed the highest overlap with evergreen type, and the lowest overlap with the deciduous group. Among the growth form types, the lowest degree of overlapping was found among the hypervolumes of trees, and arborescent-shrubs. In general, the degree of hypervolume-overlapping among different growth forms was relatively low, which indicate a trait space occupation more variable for growth forms than for environment or leaf habit.

Differences within the same groups were also remarkable (Figure 5). Only evergreen type was represented at the three environment types; it showed low functional niche space overlap (<0.3) with semideciduous (at the ridge) and (<0.2) with deciduous (at the riparian). The functional space was also very different among growth forms within the same environments (overlap < 0.4); with exception of trees and arborescent-shrubs in ridge forests, which had the highest overlap (Figure 5).

We have documented a consistent variation in functional space (based on leaf morphology and chemical composition) for woody communities distributed along the topographical gradient. Thus, communities dominated by deciduous trees with faster growth and a predominant acquisitive strategy (i.e., with high values of SLA and leaf nutrient concentration) were characteristic of riparian forests, while semideciduous shrubs and evergreen (arborescent, trees) species, characterized by a conservative strategy (i.e., with high-density leaves and a predominance of carbon-based structural defenses), dominated ridge forests. Hence, the first strong mechanism of species segregation and niche partitioning was related to habitat suitability. Along that local gradient of stress-productivity, plant establishment might be limited physiologically at one end (by stress tolerance) and by competition at the other end. Overall, leaf trait variation followed the patterns commonly associated to leaf economics spectrum (Grime, 1979; Wright et al., 2004).

This general niche differentiation began with the evolutionary accentuation of trait differences during the expansion of deciduous species on Cretaceous (Axelrod, 1966). The evolutionary and environmental drivers have determined this leaf chemical and morphological composition to improve plant functioning (Reich et al., 2003). However, we detected another trait dimension related with leaf concentration in some macronutrients (third PCA axis) that was not clearly aligned with the leaf economics spectrum. Other factors, like site differences or evolutionary history, might explain the lower values of leaf P and K associated to deciduous trees (Watanabe et al., 2007; Auger and Shipley, 2013). In a previous study on leaf chemical traits in 98 Mediterranean woody species, we also observed that leaf P and K were highly conditioned by the environment and phylogeny (de la Riva et al., accepted). These results suggest that leaf trait variation is not aligned along a single acquisition-conservation axis, and highlights the necessity of considering plant traits from other potentially independent leading dimensions such as those related with leaf chemical composition.

Within each topographical zone or environment type, the leaf biogeochemical niche partitioning would underlie the species ability to coexist by diverging on leaf nutrient composition and resource uptake. Thus, riparian communities from valley bottom showed the highest leaf biogeochemical space, which supports the existence of a greater number of functionally different species, as far as the acquisition of resources is concerned. In a regional study, the higher functional diversity (at whole plant level) was also found in wetter environments along the aridity gradient (de la Riva et al., 2017). On the contrary, the lowest diversity could be due to resource scarcity limiting the establishment of species that are not physiologically able to tolerate such abiotic constraints; as a consequence, the range of functional traits is reduced within the range of viable traits that allow plants to persist in that arid environment, in detriment of the functional space (de la Riva et al., 2016c). The more productive environments, however, are able to increase resource heterogeneity, and promote the coexistence of the species with different demand of them (de la Riva et al., 2017).

Niche partitioning within the same growth form or leaf life-span may also facilitate the coexistence among woody species. Riparian forests were dominated by deciduous tree species, which in fact was the life-span type with higher biogeochemical and morphological niche functional space. The dense shade created by the trees in this type of environment likely may induce a high among-species competition for light. However, woody species may also exhibit some degree of plasticity in nutrient uptake to respond to this competition. Plant species adapted to these productive environments, where nutrients have usually intermittent availability, might respond better to these temporal changes showing higher capacities for taking up resources and higher nutrient flexibility (Aerts, 1999; Sardans et al., 2015). Thus, our results suggest that differences in niche functional space (implying different strategies of resource capturing) among coexisting species would potentially allow to buffer competition in these resource-rich environments.

On the contrary, we detected smaller functional space in more stressful environments (i.e., ridge forests). This result could be explained because species growing in poor environments comprise traits which lead to nutrient retention and higher nutrient-use efficiency, which seem to confer a lower capacity to change their functioning in response to environmental changes (Aerts, 1999). In addition, plant size seems to be related with root expansion and nutrient uptake. Thus smaller species (more frequent in resource-poor environments) tend to show lower values of leaf nutrients than trees as a result of their lower capacity to maintain larger root systems and explore larger soil volumes to uptake nutrients (Niinemets and Kull, 2003 and references therein). As a consequence, shrub and arborescent-shrub species had lower niche space than trees. The combined effects of the slow population dynamics in poor environments (Wright, 2002), and their lower capacity to alter their elemental composition independently of the nutrient pulses (Sardans et al., 2015) would facilitate the coexistence of ecologically equivalent species.

Differences in growth form and leaf life-span may facilitate the coexistence of woody species and their niche partitioning, as showed by the lower functional overlap among them. In spite of the general patterns found along the topographical gradient (e.g., dominance of deciduous trees at riparian forests, and evergreen arborescent-shrubs at ridge forests), strong divergence were also found between their functional spaces within the same habitat type. Thus, the low degree of overlap between the functional space of deciduous trees and evergreen arborescent-shrubs at the riparian forest reflects a niche partitioning, probably promoted by light competition and resource partitioning. The dense-shaded forest understory is a relatively stressful habitat for plants (in terms of light availability), which promotes a resource gradient also related with leaf economics spectrum (Lusk and Warton, 2007). On the one hand, deciduous trees (such as Alnus glutinosa and Fraxinus angustifolia) are effective competitors, since they are able to rapidly acquire nutrients and water, that allow them to grow faster and compete more efficiently for light interception (Lambers and Poorter, 1992; Wright et al., 2001); on the other hand, evergreen shrub species (such as Rubus ulmifolius) had higher LMA, lower nutrient concentration, and are associated with shade tolerance (Williams et al., 1989; Lusk, 2002).

In the opposite extreme of the environmental gradient – at the ridge forests – woody species with semideciduous life-span seem to display more acquisitive patterns than evergreens. Shorter leaf life-span is related to higher relative nutrient requirement and lower resistance to physical hazards (Ryser, 1996; Pérez-Harguindeguy et al., 2013). In this sense, the habit of semideciduous leaf could be considered a facultative strategy, typical of shrubs in dry Mediterranean conditions (Axelrod, 1966; Zunzunegui et al., 2005). That strategy allows them to be effective competitors during favorable conditions, while they shed partly or completely their leaves during summer, so reducing water loss by transpiration (Zunzunegui et al., 2005; Ciccarelli et al., 2016; Bongers et al., 2017). Therefore, the shrub species might be segregated along the resource gradient through a trade-off between growth rate and survival. In summary, changes in plant structure such as growth form or leaf life-span, might reduce competition by differences in the ability to capture and use resources, supporting that limiting similarity was operating within the communities (Stubbs and Wilson, 2004).