Forest Ecology and Management 506 (2022) 119970 Contents lists available at ScienceDirect Forest Ecology and Management journal homepage: www.elsevier.com/locate/foreco Delayed and altered post-fire recovery pathways of Mediterranean shrubland under 20-year drought manipulation Daijun Liu a, b, c, d, Chao Zhang d, e, *, Romà Ogaya d, f, Marc Estiarte d, f, Xiwen Zhang g, Thomas A. M. Pugh b, c, h, Josep Peñuelas d, f a Department of Botany and Biodiversity Research, University of Vienna, Rennweg 14, 1030 Vienna, Austria b School of Geography, Earth and Environmental Sciences, University of Birmingham, B15 2TT Birmingham, United Kingdom c Birmingham Institute of Forest Research, University of Birmingham, B15 2TT Birmingham, United Kingdom d CSIC, Global Ecology Unit, CREAF-CEAB-CSIC, Bellaterra (Catalonia) E-08193, Spain e Optics of Photosynthesis Laboratory, Institute for Atmospheric and Earth System Research (INAR)/Forest Sciences/Viikki Plant Science Centre, Faculty of Agriculture and Forestry, University of Helsinki, PO Box 27, 00014 Helsinki, Finland f CREAF, Cerdanyola del Vallès (Catalonia) E-08193, Spain g School of Architecture and Urban Planning, Guangdong University of Technology, Guangzhou, 510090, People’s Republic of China h Department of Physical Geography and Ecosystem Science, Lund University, 22362 Lund, Sweden A R T I C L E I N F O A B S T R A C T Keywords: Increasing water deficits and severe droughts are expected to alter the dynamics of vegetation post-disturbance Succession pathways recovery by decreasing new recruitment and limiting growth in semi-arid Mediterranean ecosystems in future. Climate change However, which vegetation metrics will be shifted and how they respond over time are not clear, and the Mediterranean shrublands experimental evidence is still limited. Here we assessed the impacts of a long-term (20 years) experimental Long-term drought experiment Species diversity drought (− 30% rainfall) on the pathways of vegetation metrics related to species richness, community Community composition composition and abundance dynamics for an early-successional Mediterranean shrubland. The results indicate Abundance dynamics that the pathways of vegetation metrics were differently affected by experimental drought. The abundance of Vegetation metric space Globularia alypum follows pathway 1 (altered mature state). Simpson diversity and abundance of Erica multiflora follow pathway 2 (delayed succession) while species richness, community abundance and shrub abundance follow pathway 3 (alternative stable state). There were no significances for the resilience to extremely dry years (the ratio between the performance after and before severe events) between control and drought treatment for all vegetation metric. But, their resilience for the metrics (except Simpson diversity) to extremely dry years in 2016–17 were significantly lower than that of 2001 and of 2006–07, possibly caused by the severe water deficits in 2016–17 at mature successional stage. Principal component analysis (PCA) shows that the first two principal components explained 72.3 % of the variance in vegetation metrics. The first axis was mainly related to the changes in community abundance, shrub abundance and species richness while the second axis was related to Simpson diversity and abundance of G. alypum and E. multiflora. Principal component scores along PC1 between control and drought treatment were significantly decreased by long-term experimental drought, but the scores along PC2 were not affected. Further research should focus on successional pathways in more water-deficit conditions in Mediterranean ecosystems and the consequences of changes in vegetation recovery pathways on ecosystem functions such as biomass accumulation and soil properties. 1. Introduction having significantly altered 75% of terrestrial ecosystems, leading to abrupt biodiversity loss and ecosystem function decline (Moreno-Mateos Anthropogenic disturbances (e.g. land-use change, climate and et al., 2017; Díaz et al., 2019). Considerable efforts have been dedicated frequent fires) present an unprecedented risk for natural ecosystems, to understanding how species diversity and composition, and ecosystem * Corresponding author at: Optics of Photosynthesis Laboratory, Institute for Atmospheric and Earth System Research (INAR)/Forest Sciences/Viikki Plant Science Centre, Faculty of Agriculture and Forestry, University of Helsinki, PO Box 27, 00014 Helsinki, Finland. E-mail address: chao.x.zhang@helsinki.fi (C. Zhang). https://doi.org/10.1016/j.foreco.2021.119970 Received 13 October 2021; Received in revised form 12 December 2021; Accepted 14 December 2021 Available online 24 December 2021 0378-1127/Crown Copyright © 2021 Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). D. Liu et al. Forest Ecology and Management 506 (2022) 119970 functions recover over time (Anderson-Teixeira et al., 2013). However, have predicted that increasing aridity may reduce the post-fire recovery increasing water deficits characterized by lower soil moisture and higher (e.g. plant cover) of Mediterranean forests based on a modelling atmospheric aridity are posing greater challenges to the disturbance- approach, with large areas shifting towards open shrublands. However, recovery due to the resultant decreases in new establishment and in­ it is difficult to generalise modelling approaches on the dynamics of creases in plant mortality (Angeler and Allen, 2016; Schwalm et al., vegetation metrics for local communities, since the recovery processes 2017; Stevens-Rumann et al., 2018; Astigarraga et al., 2020). Increased greatly depend on the vegetation types, soil characteristics and biotic frequency and intensity of extreme episodes (e.g. heatwaves and severe competition (Capitanio and Carcaillet, 2008; Anderson-Teixeira et al., droughts) may trigger greater limitation on plant growth and survival, 2013; van Blerk et al., 2021b). Therefore, finding effective ways to leading to the abrupt losses in local species diversity, composition shifts capture the recovery dynamics of multiple attributes of community in and further declines in the carbon sink for terrestrial forest ecosystems response to climate change is essential. (Ciais et al., 2005; Nolan et al., 2018; Trisos et al., 2020). Therefore, Manipulative field experiments have the potential to verify the re­ accurately understanding the dynamics of ecological recovery in water- covery dynamics of species diversity, community composition and stressed ecosystems is of great importance to conserve the local diversity ecological functions and biogeochemical processes by conducting a and maintain ecosystem functions in a rapidly changing climate for continuous climate change in natural ecosystems (Wu et al., 2011; Beier coming decades. et al., 2012; Peñuelas et al., 2018). In particular, long-term (decadal) There are substantial knowledge gaps on the dynamics of ecological experiments can offer unique possibilities to detect the lagged and recovery in response to water deficit. First, which vegetation metrics accumulated impacts, which include plant physiological, demographical could be affected that is not well-addressed. Vegetation recovery is and evolutionary adjustments (Harte et al., 2015; Zhang et al., 2017; Liu mainly focused on the ecosystem functions such as vegetation biomass et al., 2020) and biotic interactions (Wright et al., 2017; Liu et al., and productivity (Anderson-Teixeira et al., 2013). Multiple vegetation 2018b). However, most of these climatic experiments have been per­ metrics such as species richness, community composition and abun­ formed in mature or equilibrium ecosystems where the structure and dance (or density) are important indicators of vegetation recovery dy­ function of plant communities are relatively complete (Beier et al., 2012; namics that are less studied (Anderson-Teixeira et al., 2013). Second, the Kröel-Dulay et al., 2015; Liu et al., 2020). Thus, vegetation responses, e. trajectories of vegetation recovery over time in a changing climate are g. species richness, community composition and productivity were not clear. Increasing summer water deficit can limit the species occur­ relatively resistant to climate treatments to the point that only minor rences of reseeding (regenerating from seed) rather than resprouting and even no net changes were detected for arctic tundra (Hudson and species, potentially shifting community composition over time (Parra Henry, 2010), temperate grassland (Grime et al., 2008) and Mediterra­ and Moreno, 2018; van Blerk et al., 2021a). However, the rates of nean shrubland (Tielbörger et al., 2014). In contrast, considerable vegetation recovery may be strengthened if some drought-adapted vegetation changes in response to climate manipulation were detected species become the dominant and reach a different mature state for the recovering communities (or early-successional status) since the (Anderson-Teixeira et al., 2013). In addition, severe droughts may delay species assembles were unstable (Prieto et al., 2009; Liu et al., 2017; or shift vegetation post-recovery via the alteration in growth resilience Parra and Moreno 2018; van Blerk et al., 2021a). It fits well with hier­ (the capacity to reach pre-disturbance performance) (Parra and Moreno, archical theory that abrupt and nonlinear changes in community 2018; Lloret et al., 2011, 2012). Third, effective approaches to investi­ structure would likely occur over the long term by species reordering gate ecological recovery at stand level is greatly limited. Ecosystem and competition outcomes (Smith et al., 2009; Luo et al., 2011). In modelling has been applied to estimate the recovery time of global addition, the impact of extreme natural droughts on community and forests after disturbances such as biomass harvest, large fires and insect population dynamics could be further intensified, due to changes in the outbreaks (Pugh et al., 2019), but the effect of increasing drought stress resilience in the climate-manipulated ecosystems (Barbeta et al., 2015; on recovery was not investigated. Interestingly, Baudena et al. (2020) Liu et al., 2015, 2018a; Parra and Moreno 2018). But it has also been Fig. 1. Conceptual diagram hypothesis­ ing four successional pathways in which climate change may impact the dy­ namics of vegetation recovery. The blue line represents the value of a vegetation metric following fire disturbance in the absence of experimental drought (con­ trol). Metrics include species richness, community composition and abundance dynamics (for community, shrub and species levels). The yellow lines indicate the four pathways simulated by experi­ mental drought (− 30% rainfall), repre­ senting a possible impact of climate change (Fig. 1a). The vertical dashed line separated the vegetation recovery into two stages (young and mature) in control, the peak year was selected based on the highest value of vegetation metrics. The right figures (Fig. 1b,c,d,e) show the percent changes of vegetation metrics (calculated as: 100*(drought - control)/control), with the horizontal lines indicate the values at the beginning of experiment. Therefore, long-term drought experiment and extreme droughts could lead to four different succession pathways for the early successional ecosystem. 1) Altered mature state (pathway 1): the species can obtain resource benefits under drought conditions and perform better than the control (Fig. 1b). 2) Delayed succession (pathway 2): the new estab­ lishment and growth are transiently limited, and the recovery of vegetation metrics may take a longer time to reach its expected mature state (Fig. 1c). 3) Alternative stable state (pathway 3): some species are permanently filtered out due to species reordering, plant community will achieve an alternative stable status (Fig. 1d). 4) Advanced degradation (pathway 4): species establishment and growth are strongly limited and tend to degrade over time due to the accumulated and legacy effects of drought (Fig. 1e). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 2 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 reported that vegetation may have the resilience after extreme droughts period of 1951–2018 was 571.2 mm, with nearly 62 % occurring in via the stabilisation process of increasing new recruitment and regrowth spring and autumn (354.1 mm). Mean annual temperature (MAT) dur­ (Seifan et al., 2010; Lloret et al., 2012; Ploughe et al., 2019). To date, ing the same period was 15.3 ◦ C, and the hottest season occurred in few long-term drought experiments have been performed in recovering summer (mean temperature is 22.8 ◦ C) (Fig. S1). The soil type is and unstable habitats; detailed knowledge on the succession pathways calcareous and composed of marlstone and limestone. The soil depth of species diversity, composition and abundance dynamics, is highly above the bedrock is about 10–40 cm. The study site has suffered two limited. severe wildfires, one in summer 1982 and the other in spring 1994. The The successional recovery is strongly affected by drought in the conifer forest (dominated by Pinus halepensis) was damaged by the two water-limited Mediterranean-type ecosystems. The Mediterranean Basin fires. Thus, most of the current vegetation, mainly Mediterranean is one of the most highly disturbed regions, due to a combination of shrubs, has sprouted from the underground organs after the fires. The anthropogenic and natural activities, and it has been demonstrated that vegetation in this site has the potential to recover to become a secondary these disturbed (or successional) communities are sensitive to ongoing forest after the two fires. But tree seedlings (e.g. Pinus halepensis) were and future climate change (especially drought) (Myers et al., 2000; cleared irregularly to maintain the vegetation as shrubland. The height Doblas-Miranda et al., 2015; Cramer et al., 2018). Negative impacts on of shrubland vegetation is<1.5 m. The canopy vegetation cover is about ecological recovery (e.g. forest degradation) are expected to be further 40–60%. There are nearly 30 Mediterranean shrub and herbaceous exacerbated by the projected 2–4 ◦ C summer temperature increase and species in this study site. The dominant shrubs are Globularia alypum continuous 30% summer precipitation decrease by the end of this cen­ (GA) and Erica multiflora (EM). The abundance of GA and EM relative to tury (Cramer et al., 2018; Lionello and Scarascia, 2018). Experimental total abundance is larger than 50% in our study site. Moreover, there are results have indicated that drought has already delayed the post-fire other shrub species, such as Rhthtosmarinus officinalis (RO), Ulex parvi­ successional process by decreasing the species richness and density of florus (UP) and Dorycnium pentaphyllum (DP). The undergrowth also recruited seedlings (Lloret et al., 2004, 2009; del Cacho et al., 2013). includes some herbaceous species (mostly Poaceae species (Poaceae Once the communities reached a more stable status, increases in species spp.)). richness under drought can happen, due to the establishment of drought- adapted species (e.g. G. alypum) (Liu et al., 2018b). But abundance 2.2. Field manipulations changes may vary greatly with functional groups (e.g. herbs vs shrubs). For instance, in an early-successional Mediterranean shrubland ecosys­ The treatments of control and drought experiments were initiated in tems, both long-term warming and drought had no significant influence 1999 (1998 was a pre-treatment year). Three blocks with two plots (4 m on the abundance of shrubs (Liu et al., 2017). Indeed, the shrubs and the × 5 m) each were selected along a south-facing slope (ca. 13%). Treat­ dominant species (G. alypum and E. multiflora) may acquire resource ments were randomly distributed in each block (Beier et al., 2004; Liu benefits due to their long life-history, physiological adaptions (e.g. et al., 2017, 2020). Scaffolding (1.2 m height) was constructed in all of increased photosynthesis and water-use efficiency) (Liu et al., 2016; the plots. But there were no curtains or covers in control plots. Future Zhang et al., 2017) and delayed leaf-unfolding and flowering (Peñuelas drought condition was simulated by intercepting and removing rainfall et al., 2004; Prieto et al., 2008). However, the successional pathways of with transparent waterproof covers during the wet seasons (spring and species diversity, community composition and abundance dynamics autumn). The drought curtains were activated by a rain sensor and affected by long-term drought, are still unclear. covered the drought plots during rainfall events larger than 5 mm and We conducted a drought (− 30% rainfall) manipulation in an early- were retracted once the rain ended. The automatic system for drought successional Mediterranean shrubland to study the recovery dynamics was also damaged by lighting in 2015 and then was redesigned to pas­ over the period of 1999–2018. The study site had suffered a wildfire 4 sive continuous rain interception using transparent V-type plastic over years before the drought manipulation. The amount of rainfall reduction the scaffolds designed to intercept nearly 30% of the annual rainfall in our experiment fell into the range of rainfall changes in the projection (Fig. S2), the similar percentage that was intercepted with the automatic of climatic models (up to − 30% with the scenario of 2–4 ◦ C) for the system. The cover system had minor influence on the soil and surface Mediterranean Basin during the coming decades (Cramer et al., 2018; soil temperature. Soil moisture at 0–15 cm depth was measured Lionello and Scarascia, 2018). We hypothesized that long-term (20- biweekly using Time Domain Reflectometry (TDR, Tektronix 1502C; year) climatic manipulation, combined with severe droughts (a 10th Tektronix, Beaverton, OR, USA) in each plot. Throughout the study percentile threshold of annual precipitation), would shift the pathways period of 1999–2018, the 30% rainfall reduction in drought plots of vegetation recovery and reduce the vegetation resilience to severe reduced the soil moisture at 0–15 cm depth by 16.7%, compared to events in this Mediterranean early-successional shrubland. The objec­ control plots (15.50 ± 0.7% and 18.59 ± 0.6% for drought and control tives were to 1) study the impacts of long-term drought on species di­ respectively, n = 20). versity, community composition, abundance dynamics (for community, shrub and two dominant species; respectively) and to identify the 2.3. Vegetation measurements pathways (altered mature state, delayed succession, alternative stable state and advanced degradation) explained in Fig. 1; 2) compare the Vegetation was measured by a nondestructive pin-pointing method likely differences in the resilience to severe droughts between drought conducted annually during the non-growing season in late July and the treatment and control; 3) analyze the effect of experimental drought on beginning of August from 1998 (pre-treatment) to 2018. With the pin- the changes in spatial variation of the vegetation metrics, which can be point method, we successfully recorded every hit within the vegeta­ investigated by the change in the main principal components. tion profile that allowed us to estimate total number and density fre­ quencies for every species and that were used to follow the changes over 2. Materials and methods time in the plant community structure and composition. We conducted vegetation measurements since the pre-treatment year in five parallel 3- 2.1. Study site m transects (0.8 m apart) marked permanently in the central area of each plot (3 m × 4 m). Each transect had a total of 61 sampling points The study site is located in the Garraf Natural Park south of Barce­ that were separated every 5 cm, totaling 305 points per plot. For the lona, northeastern Spain (41◦ 18́ N, 1◦ 49́ E; 210 m a.m.s.l.). The climate is sampling, we lowered a 2 m long steel pin (3 mm diameter) through typically Mediterranean-type, characterised with dry and hot summers each sampling point and recorded, for each plant hit, the organ (e.g. leaf, (June to August) and wet springs (March to May) and autumns branch and reproductive of flower and fruit) and its status (alive or (September to November). Mean annual precipitation (MAP) during the dead). Our measurements were performed during the non-growing 3 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Fig. 2. The changes in (a) species richness and (b) composition diversity in control and drought treatments throughout the study period of 1998–2018. The year of 1998 is the pre-treatment year and the years with red arrows indicated the extreme dry years during the study period. Vertical bars represent the standard errors of the means (n = 3). Shaded grey bars are the extremely dry years of 2001, 2006–07 and 2016–17. The solid white dots and dashed lines at the bottom (c, d) indicate the percentage difference between control and drought plots (calculated as: 100* (drought - control)/control). The dashed grey horizontal lines (y = 0) separate the difference between positive and negative values. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) season, the dead hits included both those resulting from naturally Climate variables of monthly mean temperature and total precipi­ drought and those caused by the climate manipulation. Since the dead tation for 1951–2018 were compiled from a nearby meteorological hits occurred in both control and drought plots. station during the period of 1951–2010 and our study site during the study period of 1998–2018 (linear regression for the two sites for 2.4. Vegetation metrics and climate variables 1998–2010 of R2 = 0.99 for monthly temperature and R2 = 0.83 for monthly precipitation). Based on this profile, we calculated the histor­ Species richness and abundance dynamics were measured yearly in ical climate variables of annual mean temperature and total precipita­ each plot during the study period of 1998–2018. The abundance (hits) of tion for the study site. The natural extremely dry years were obtained all species in the plot was recorded and used to calculate the index of according to the lowest 10th percentile precipitation year during the abundance dynamics and Simpson diversity. Here abundance metric period of 1951–2018 (Knapp et al., 2015). During the experimental was calculated as the number of alive hits or frequency. The number of period of 1998–2018, extremely dry years occurred in 2001, 2006–07 alive hits for each species were naturally log-transformed and then and 2016–17 (Fig. S1). summed to get the equal weight of community abundance. The shrub abundance each plot was obtained by summing all the log-transformed 2.5. Statistical analyses hits of shrub species according to the growth-form identified in our study site. The abundance of the two dominant species G. alypum and We selected generalized additive mixed models (GAMMs) to analyze E. multiflora was analyzed separately. Community composition was the dynamics of vegetation metrics with a smoothing function using the estimated by the Simpson diversity index for each plot throughout the mgcv package (version 1.8.3) (Wood, 2017), because vegetation metrics study period. It was calculated as: of post-fire recovery process are nonlinear and flexible. The bam func­ ∑ n(n − 1) tion was implemented in the GAMM analysis. The smoothing functions D = 1− N(N − 1) could be used for testing potentially different trends over time for different treatments (control, drought). The model formula was per­ where n is the number of alive hits in a specific species and N is total formed with ordered factors that included intercept difference and number of alive hits for all species. Simpson diversity takes into account reference smoothing: vegetation metric ~ treatment + s(year) + s(year, the number of plant species and the abundance of each species in a by = treatments). Here, treatment that had two levels, drought and community. The values of its range from 0 to 1, with higher values control, was considered as a reordered factor, and year as a numeric indicating higher diversity. Simpson diversity were calculated with the variable from 1998 (pre-treatment year) to 2018. The plot nested within vegdist function in the vegan 2.5–6 package (Oksanen et al., 2017) in R block was treated as the random factor in the analyses. In order to version 3.3.0. explain the temporal autocorrelation in the time series, we specified the 4 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Fig. 3. Changes in abundance at community and shrub levels (a, b) throughout the study period of 1998–2018. Abundances were calculated as the natural log of hits for all the species in the plot and all shrub species, respectively. The year of 1998 is the pre-treatment year and the years with red arrows indicated the extreme dry years during the study period. Vertical bars represent the standard errors of the means (n = 3). Shaded grey bars are the extremely dry years of 2001, 2006–07 and 2016–17. Differences are shown in the bottom (c, d) as in Fig. 2. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) correlation structure by adding year in the model analysis. The model study years were analyzed using the Performance Analytics package result can imply the differences between treatments and the two smooth (version 2.0.4) in R (Peterson et al., 2019). Here, we used the packages terms: reference smooth for the control and difference smooth (control- FactoMineR (version 1.34) for the analysis (Husson et al., 2013) and drought) over time (the model result for all vegetation metrics is factoextra (version 1.0.7) for ggplot2-based visualization (Kassambara, detailed in Table S1). 2016). The scores of the vegetation variables in each plot during the We conducted the assessments on the resilience to extremely dry period 1999–2018 were obtained according to the score function. Here years of 2001, 2006–07 and 2016–17. Resilience is the capacity to reach we only selected the scores along the first two dimensions. Due to the pre-disturbance performance, which was calculated as the values after non-independence of the measurements in control and drought over an extremely dry year versus those in the pre-drought year in each plot time, we included an autocorrelation structure in the model using nlme (PostDry / PreDry) (Lloret et al., 2011). There were two continuous package. The significant differences for the interaction of treatments extremely dry years in 2006–07 and 2016–17. The PreDry years are (control and drought) and year were analysed. Plot was selected as the 2000, 2005 and 2015 while PostDry years are 2002, 2008 and 2018. The random factor in the best model. Repeated measure for a two-way three extremely dry events occurred at different stages of recovery: analysis of variance (ANOVA) was applied to test the significance. The young and mature. These two stages were separated by the peak year of mean score of each treatment × year combination was calculated by the highest value of community abundance (2010). Two-way ANOVA groupwiseMean function from the rcompanion package, and the confi­ (treatment and extreme dry years) was used to test the significant in­ dence interval for the mean scores which were presented in the Fig. 6b, fluences. Once the significance was obtained, we further analysed the c. ANOVA with Tukey’s HSD (honest significant difference) post-hoc test. All the analyses were performed in R (version 4.0.5). Principal component analysis (PCA) was conducted to explore and visualize the relationships between the multiple quantitative vegetation 3. Results metrics and to explore which metrics are the most crucial in dis­ tinguishing the functional differences for treatments. The vegetation 3. 1 Species richness and Simpson diversity in response to long-term metrics of species richness, Simpson diversity, community abundance, drought shrub abundance, abundance of G. alypum (GA) and E. multiflora (EM) Species richness and Simpson diversity were significantly decreased were selected. We applied Horn’s parallel analysis (paran package of by drought (both p < 0.001) during the period 1999–2018 based on the version 1.5.2) to determine how many axes of vegetation metrics should GAMM test (Fig. 2a and Fig. 2b; respectively) (Table S1). The smoothing be retained in principal component analysis. The result of the parallel in control for species richness and Simpson diversity were nonlinear, analysis indicated that all of the 6 metrics should be retained in the PCA and they were significantly different from 0 in the model analysis. We do (Fig. S3). The treatment (control and drought) was considered as a not observe a significant interaction of treatment and year for these two categorical variable. The relationship between the vegetation metric and vegetation metrics. Based on the result of model analysis and percentage 5 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Fig. 4. The response of abundance (a, b) of two dominant species (G. alypum and E. multiflora) throughout the study period of 1998–2018. The abundance was calculated as the natural log of hits for the species in the plot. The year of 1998 is the pre-treatment year. Vertical bars represent the standard errors of the means (n = 3). Shaded grey bars are the extremely dry years of 2001, 2006–07 and 2016–17. Differences are shown in the bottom as in Fig. 2. difference between control and drought, we can conclude that species 3.3. The variation in vegetation metric space over time richness and Simpson diversity follow pathway 3 (alternative stable state) and pathway 2 (delayed succession), respectively. The first two principal components (PC; axes) of the PCA explained 72.3 % of the variance in vegetation metrics, PC1 (47.4 %) and PC2 (24.9 %) (Fig. 6a). The first dimension (PC1) was related to the changes 3.1. Abundance dynamics affected by drought treatment in community abundance, shrub abundance and species richness, that were highly positively associated in the first dimension (contributed by Both community abundance and shrub abundance were significantly 56.3 %, 52.9 % and 48.1 %, respectively). However, the second decreased by drought treatment over the study period (Fig. 3a,b; dimension (PC2) was linked to the Simpson diversity (contributed by Table S1). The GAMM tests suggest that the smoothing was significantly 47.3 %) and abundance of G. alypum and E. multiflora (contributed by increasing over the study period (both p < 0.001). The changes in 74.2 % and a 37.2 %, respectively). Simpson diversity was negatively abundance at community and shrub levels by drought follow the correlated with the changes in the abundance of G. alypum and pathway 3 (alternative stable state) (Fig. 3c, d). E. multiflora (Fig. S4). In addition, the scores for the two PCs varied Drought treatment had different effects on the abundance for the two differently with the treatments over time (Fig. 6b). PC1 scores were dominant species (Fig. 4a, b; Table S1). Drought treatment significantly significantly decreased by experimental drought (anova, p < 0.001). increased abundance of G. alypum (p < 0.001) and there was an inter­ PC2 scores were not affected by treatment, but the changes varied active effect of drought and year (p < 0.05). Yet it did not affect abun­ significantly over time (p < 0.001). dance of E. multiflora during the study period. The changes in abundance of G. alypum and of E. multiflora support pathway 1 (altered mature 4. Discussion state) and pathway 2 (delayed succession); respectively. 4.1. Changes in the species richness and community diversity 3.2. The resilience to extreme droughts Long-term drought had significant impacts on species richness and composition. Species richness and Simpson diversity were significantly We observed that there were no significant differences for the resil­ decreased by drought treatment but the smoothing differences between ience of any of the vegetation metrics between control and drought control and drought did not change over time. The changes in species treatments based on the ANOVA test (Fig. 5; Table S2). However, there richness under drought fit well with the pathway 3 (altered maturity were significant differences among the three extreme events for the state) that reach a relatively stable state over time (Fig. 2a), possibly due resilience of vegetation metrics (except Simpson diversity). Moreover, to the loss of drought-sensitive species (e.g. Galium lucidum, Ulex par­ the values of resilience were similar between the first two extreme viflorus and Lithodora fruticosa). Simpson diversity in drought appears to events (2001 and 2006–07), which occurred during the early succes­ follow pathway 2 (delayed succession) that took a longer time to reach sional stage, and the values at the mature successional stage (2016–17) its expected mature state. This can be attributed to increases in were lower than the values at the early stage. 6 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Fig. 5. The resilience (Rs) in response to extremely dry years. The values of species richness (a), Simpson diversity (b), community abundance (c), shrub abundance (d), abundance of G. alypum (e) E. multiflora (f) are shown. The extremely dry years were 2001, 2006–07 and 2016–17; respectively. Blue and yellow indicate the control and drought plots. The different letters indicate the significance among the extremely dry years by the ANOVA with Tukey’s HSD post-hoc test. And no letters in the plots mean no significance was found in the analysis. Vertical bars represent the standard errors. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) abundance of drought-resistant species (e.g. G. alypum, Polygala rupestris recover as the mature state in control plots and reach a relatively stable and Fumana thymifolia) in dry condition (Liu et al., 2018b). However, we state. Less abundance (density) for both community and shrub levels observe the resilience for species richness to the extreme dry years of under drought might be associated with the limitation on occurrence of 2016–17 (mature stage) was lower than to those of 2001(young stage). new organs (e.g. leaves, branch and flowers) (Prieto et al., 2009; Liu This result suggests more negative influences by extreme events during et al., 2017). However, the significant decrease in both community the mature stage of succession, although this response could also be abundance and shrub abundance over time, especially in the extremely related to the cumulative effect of multiple droughts on the recovering dry years of 2016–17. One reason may be related to self-thinning once vegetation. This contrasts with the response reported for mature Medi­ the plant communities were approaching mature status (Yu et al., 2019). terranean shrubland in Israel where no net changes in species richness The resilience of the community abundance and shrub abundance to and biomass accumulation were reported under decadal drought (− 30% extreme events were decreased over time; the values were lower in rainfall), partly due to high resistance of mature communities mature stage (in 2016–17) than in young stage (in 2001 and 2006–07) (Tielbörger et al., 2014). Moreover, it was reported that few sites had (Fig. 5c, d). Another reason that the rainfall of these two years was significant changes in shrubland species richness by experimental relatively low (average 402 mm), which was much less than the average drought over time for study sites across the climate gradients of Europe rainfall (586 mm across 1951–2018). (Kröel-Dulay et al., 2015). The drought impacts on abundance for the dominant species differed greatly. Drought significantly increased the abundance of G. alypum compared to control, and this drought difference was increased over 4.2. Experimental drought affected abundance dynamics over time time. But experimental drought decreased abundance of E. multiflora in the first few years and then had almost no influences after 2007. Community abundance and shrub abundance were significantly Therefore, the abundance of G. alypum follows pathway 1 (alternative decreased by experimental drought. These two metrics under drought fit stable state), but E. multiflora follows the pathway 2 (delayed with the pathway 3 (alternative stable state) that being restricted to 7 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Fig. 6. The variation of vegetation metric space (a) as explained by the first two dimensions (b, c) of PCAs. The blue dots indicate control and yellow dots indicate drought plots. The changes for PC1 score were analysed by the interaction of treatment (control and drought) and year using the gls function in nlme package. Plot was selected as the random factor in the best model. The mean score of each treatment × year combination was calculated by groupwiseMean function from rcompanion package. The blue lines show the control, while orange lines show drought. The error bars represented the confidence interval of the plots. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) succession). Previous studies have reported that the dominant species manipulations to understand the real types of succession pathway over (especially G. alypum) have advantages to endure warmer and drier time. For pathway 4 (advanced degradation) to occur, the drought conditions by adjusting their photosynthesis and water-use efficiency resistance and plasticity are likely overwhelmed (Tielbörger et al., (Liu et al., 2016; Zhang et al., 2017) and phenological alteration (e.g. 2014), which may take many years or decades to become apparent. delaying or advancing the leaf appearance) (Peñuelas et al., 2004; Prieto Thus, it is plausible that given a longer drought perturbation or et al., 2008). However, the abundance of G. alypum abruptly declined in moderately more severe extremes, these systems may yet shift towards both control and drought plots in extremely dry years of 2016–17, which pathway 4, with more severe consequences on ecosystem functions in may be caused by accumulatively water deficits. We also observed a future. Thus, continued long-term drought manipulation studies are lower resilience in the period of 2016–16 than that of extremely dry essential to provide an early warning of how climate change affects years of 2001 and 2006–07 (Fig. 5g, h). Thus, more severe extreme species diversity, community composition and abundance dynamics for drought in future may decrease the growth resilience of these two spe­ recovering communities. cies. Thus, the two dominant shrub species may lose competitive ad­ These impacts on recovery also need to be put in the context of vantages in hot and/or dry conditions if other more drought-adapted changing disturbance regimes. Global warming and intensified water species (e.g. Polygala rupestris and Fumana thymifolia). stress may increase the occurrence and severity of fires in the Mediter­ ranean ecosystems (Díaz-delgado et al., 2002; Doblas-Miranda et al., 4.3. Change in vegetation metric space over time 2015; Cramer et al., 2018). This would shift larger areas of the Medi­ terranean Basin into an earlier successional stage, exposing them to the The first two principal components clearly reflect the variation of the altered successional processes that we have described herein, and thus vegetation metric space over time. In the analysis, the first axis was accelerating a shift in species diversity and community composition largely driven by the changes in community abundance, shrub abun­ across the region. Changes in disturbance types and intensity may also dance and species richness (Fig. 6a). High positive correlations were be important. Whilst Mediterranean-type ecosystems typically have a observed for these three metrics in the positive direction of first fast recovery after fire disturbance (Pausas et al., 2008; Kelly et al., dimension (Fig. S4). The second axis was mainly driven by the Simpson 2020), many studies have indicated that high burn severity appears diversity and abundance of G. alypum and E. multiflora, and Simpson likely to decrease the ability of Mediterranean forest to recover (Fer­ diversity was negatively correlated with the changes in abundance of nandez-Manso et al., 2016). Increases in fire severity may interact with two dominant species (Fig. 6a). However, the scores along PC1 were the selection pressure of drought to further alter Mediterranean eco­ significantly higher in control than drought, indicating the vegetation systems. Moreover, habitat conditions such as soil water availability, space was shifted by experimental drought over time. The assessment of nutrient content and soil seed bank as well as occurrence of alien species variation in vegetation metrics could be useful to understand the metric may also affect the post-disturbance of local species diversity and community composition (Doblas-Miranda et al., 2015; van Blerk et al., space associated with plant ecological strategies (Díaz et al., 2016; Bruelheide et al., 2018). 2021a,b; Keeley and Brennan, 2012). Assessing the combined effects of these multiple dimensions of change is a major challenge, but our results highlight the important role of drought stress during the recovery pro­ 4.4. Greater structure and composition shift for future climate scenarios cess in leading to long-term alterations in species diversity, composition and abundance dynamics. Through their destabilising role within eco­ The changes of vegetation metrics in drought condition varied systems, disturbances bring plant communities to a state in which they greatly over time, it is still challenging to definitively attribute to are most sensitive to climate change. pathway types. Such attribution should become easier as the timeseries is extended, highlighting the importance of conducting multi-decadal 8 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Declaration of Competing Interest risks to sustainable development in the Mediterranean. Nat. Clim. Change 8, 972–980. Díaz-delgado, A.R., Lloret, F., Pons, X., Terradas, J., 2002. Satellite evidence of The authors declare that they have no known competing financial decreasing resilience in Mediterranean plant communities after recurrent wildfires. interests or personal relationships that could have appeared to influence Ecology 83, 2293–2303. the work reported in this paper. Díaz, S., Kattge, J., Cornelissen, J.H.C., Wright, I.J., Lavorel, S., Dray, S., Reu, B., Kleyer, M., Wirth, C., Prentice, I.C., et al., 2016. The global spectrum of plant form and function. Nature 529, 167–171. Díaz, S., Settele, J., Brondízio, E., 2019. Report of the Plenary of the Intergovernmental Acknowledgements Science-Policy Platform on Biodiversity and Ecosystem Services on the work of its seventh session. Intergovernmental Science-Policy Platform on Biodiversity and We greatly thank all the people who maintained the field climate Ecosystem Services (IPBES) 7, 45. experiments and conducted the samplings. Funding for this research was Doblas-Miranda, E., Martínez-Vilalta, J., Lloret, F., Álvarez, A., Ávila, A., Bonet, F.J., Brotons, L., Castro, J., Curiel Yuste, J., Díaz, M., et al., 2015. Reassessing global provided by the Spanish government project PID2019-110521GB-I00, change research priorities in mediterranean terrestrial ecosystems: How far have we Catalan government project SGR2017-1005, Academy of Finland come and where do we go from here? Glob. Ecol. Biogeogr. 24, 25–43. (340744, SensBVOCs), Fundación Ramón Areces grant ELEMENTAL- Fernandez-Manso, A., Quintano, C., Roberts, D.A., 2016. Burn severity influence on post- fire vegetation cover resilience from Landsat MESMA fraction images time series in CLIMATE, the European Research Council (ERC) under the European Mediterranean forest ecosystems. Remote Sens. Environ. 184, 112–123. Union’s Horizon 2020 research and innovation program (grant agree­ Grime, J.P., Fridley, J.D., Askew, A.P., Thompson, K., Hodgson, J., Bennett, C., 2008. ment No. 758873, TreeMort), and the FWF Austria science Fund (Lise Long-term resistance to simulated climate change in an infertile grassland. PNAS 105, 10028–10032. Meitner Programme M2714-B29). Harte, J., Saleska, S.R., Levy, C., 2015. Convergent ecosystem responses to 23-year ambient and manipulated warming link advancing snowmelt and shrub encroachment to transient and long-term climate-soil carbon feedback. Glob. Change Author contributions Biol. 21, 2349–2356. Hudson, J.M.G., Henry, G.H.R., 2010. High Arctic plant community resists 15 years of D.L. and C.Z. developed the ideas and analysed the data; J.P. and R. experimental warming. J. Ecol. 98, 1035–1041. Husson, F., Josse, J., Le, S., Maintainer, J.M., 2013. Multivariate exploratory data O. conceived the ideas and collected the data; D.L and C.Z wrote the analysis and data mining with R. R package version, 1.1.29. manuscript with contributions from all co-authors. Kassambara, A., 2016. Practical guide to principal component methods in R: PCA, M (CA), FAMD, MFA, HCPC, factoextra. Keeley JE & Brennan TJ. 2012. Fire-driven alien invasion in a fire-adapted ecosystem. Data availability Oecologia 169: 1043-1052. Kelly, L.T., Giljohann, K.M., Duane, A., Aquilué, N., Archibald, S., Batllori, E., Bennett, A. F., Buckland, S.T., Canelles, Q., Clarke, M.F., et al., 2020. Fire and biodiversity in the All data and code are available upon request to CZ (chao.x.zhan­ Anthropocene. Science 370, 6519. g@helsinki.fi). Knapp, A.K., Hoover, D.L., Wilcox, K.R., Avolio, M.L., Koerner, S.E., La Pierre, K.J., Loik, M.E., Luo, Y., Sala, O.E., Smith, M.D., 2015. Characterizing differences in precipitation regimes of extreme wet and dry years: Implications for climate change Appendix A. Supplementary material experiments. Glob. Change Biol. 21, 2624–2633. Kröel-Dulay, G., Ransijn, J., Schmidt, I.K., Beier, C., De Angelis, P., de Dato, G., Dukes, J. Supplementary data to this article can be found online at https://doi. S., Emmett, B., Estiarte, M., Garadnai, J., et al., 2015. Increased sensitivity to climate change in disturbed ecosystems. Nat. Commun. 6, 1–7. org/10.1016/j.foreco.2021.119970. Lionello, P., Scarascia, L., 2018. The relation between climate change in the Mediterranean region and global warming. Reg. Environ. Change 18, 1481–1493. References Liu, D., Estiarte, M., Ogaya, R., Yang, X., Peñuelas, J., 2017. Shift in community structure in an early-successional Mediterranean shrubland driven by long-term experimental warming and drought and natural extreme droughts. Glob. Change Biol. 23, Anderson-Teixeira, K.J., Miller, A.D., Mohan, J.E., Hudiburg, T.W., Duval, B.D., 4267–4279. DeLucia, E.H., 2013. Altered dynamics of forest recovery under a changing climate. Liu, D., Llusia, J., Ogaya, R., Estiarte, M., Llorens, L., Yang, X., Peñuelas, J., 2016. Glob. Change Biol. 19, 2001–2021. Physiological adjustments of a Mediterranean shrub to long-term experimental Angeler, D.G., Allen, C.R., 2016. Quantifying resilience. J. Appl. Ecol. 53, 617–624. warming and drought treatments. Plant Sci. 252, 53–61. Astigarraga, J., Andivia, E., Zavala, M.A., Gazol, A., Cruz-Alonso, V., Vicente-Serrano, S. Liu, D., Ogaya, R., Barbeta, A., Yang, X., Peñuelas, J., 2015. Contrasting impacts of M., Ruiz-Benito, P., 2020. Evidence of non-stationary relationships between climate continuous moderate drought and episodic severe droughts on the aboveground- and forest responses: Increased sensitivity to climate change in Iberian forests. Glob. biomass increment and litterfall of three coexisting Mediterranean woody species. Change Biol. 26, 5063–5076. Glob. Change Biol. 21, 4196–4209. Barbeta, A., Mejía-Chang, M., Ogaya, R., Voltas, J., Dawson, T.E., Peñuelas, J., 2015. The Liu, D., Ogaya, R., Barbeta, A., Yang, X., Peñuelas, J., 2018a. Long-term experimental combined effects of a long-term experimental drought and an extreme drought on drought combined with natural extremes accelerate vegetation shift in a the use of plant-water sources in a Mediterranean forest. Glob. Change Biol. 21, Mediterranean holm oak forest. Environ. Exp. Bot. 151, 1–11. 1213–1225. Liu, D., Peñuelas, J., Ogaya, R., Estiarte, M., Tielbörger, K., Slowik, F., Yang, X., Baudena, M., Santana, V.M., Baeza, M.J., Bautista, S., Eppinga, M.B., Hemerik, L., Bilton, M.C., 2018b. Species selection under long-term experimental warming and Garcia, M.A., Rodriguez, F., Valdecantos, A., Vallejo, V.R., et al., 2020. Increased drought explained by climatic distributions. New Phytol. 217, 1494–1506. aridity drives post-fire recovery of Mediterranean forests towards open shrublands. Liu, D., Zhang, C., Ogaya, R., Estiarte, M., Peñuelas, J., 2020. Effects of decadal New Phytol. 225, 1500–1515. experimental drought and climate extremes on vegetation growth in Mediterranean Beier, C., Beierkuhnlein, C., Wohlgemuth, T., Peñuuelas, J., Emmett, B., Korner, C., de forests and shrublands. J. Veg. Sci. 31, 1–12. Boeck, H., Christensen, J.H., Leuzinger, S., Janssens, I.A., et al., 2012. Precipitation Lloret, F., Escudero, A., Iriondo, J.M., Martínez-Vilalta, J., Valladares, F., 2012. Extreme manipulation experiments-chanllenges and recommendations for the future. Ecol. climatic events and vegetation: The role of stabilizing processes. Glob. Change Biol. Lett. 15, 899–911. 18, 797–805. Beier, C., Emmett, B., Gundersen, P., Tietema, A., Estiarte, M., Gordon, C., Gorissen, A., Lloret, F., Keeling, E.G., Sala, A., 2011. Components of tree resilience: effects of Llorens, L., Roda, F., Williams, D., 2004. Novel approaches to study climate change successive low-growth episodes in old ponderosa pine forests. Oikos 120, effects on terrestrial ecosystems in the field : drought and passive nighttime 1909–1920. warming. Ecosystems 7, 583–597. Lloret, F., Peñuelas, J., Prieto, P., Llorens, L., Estiarte, M., 2009. Plant community Bruelheide, H., Dengler, J., Purschke, O., Lenoir, J., Jiménez-Alfaro, B., Hennekens, S.M., changes induced by experimental climate change: Seedling and adult species Botta-Dukát, Z., Chytrý, M., Field, R., Jansen, F., et al., 2018. Global composition. Perspect. Plant Ecol. Evol. Systemat. 11, 53–63. trait–environment relationships of plant communities. Nat. Ecol. Evol. 2, Lloret, F., Siscart, D., Dalmases, C., 2004. Canopy recovery after drought dieback in 1906–1917. holm-oak Mediterranean forests of Catalonia (NE Spain). Glob. Change Biol. 10, del Cacho, M., Peñuelas, J., Lloret, F., 2013. Reproductive output in Mediterranean 2092–2099. shrubs under climate change experimentally induced by drought and warming. Luo, Y., Melillo, J., Niu, S., Beier, C., Clarks, J.S., Classen, A.T., 2011. Coordinated Perspect. Plant Ecol. Evol. Systemat. 15, 319–327. approaches to quantify long-term ecosystem dynamics in response to global change. Capitanio, R., Carcaillet, C., 2008. Post-fire Mediterranean vegetation dynamics and Glob. Change Biol. 17, 843–854. diversity: A discussion of succession models. For. Ecol. Manage. 255, 431–439. Moreno-Mateos, D., Barbier, E.B., Jones, P.C., Jones, H.P., Aronson, J., López-López, J.A., Ciais, P., Reichstein, M., Viovy, N., Granier, A., Ogée, J., Allard, V., Aubinet, M., McCrackin, M.L., Meli, P., Montoya, D., Rey Benayas, J.M., 2017. Anthropogenic Buchmann, N., Bernhofer, C., Carrara, A., et al., 2005. Europe-wide reduction in ecosystem disturbance and the recovery debt. Nat. Commun. 8, 8–13. primary productivity caused by the heat and drought in 2003. Nature 437, 529–533. Myers, N., Mittermeier, R.A., Mittermeier, C.G., da Fronseca, G.A.B., Kent, J., 2000. Cramer, W., Guiot, J., Fader, M., Garrabou, J., Gattuso, J.P., Iglesias, A., Lange, M.A., Biodiversity hotspots for conservation priorities. Nature 403, 853–858. Lionello, P., Llasat, M.C., Paz, S., et al., 2018. Climate change and interconnected 9 D. Liu et al. Forest Ecology and Management 506 (2022) 119970 Nolan, C., Overpeck, J., Allen, J.R.M., Anderson, P., Betancourt, J.L., Binney, H.A., Seifan, M., Tielbörger, K., Kadmon, R., 2010. Direct and indirect interactions among Simon, B., Bush, M.B., Chase, B.M., 2018. Past and future global transformation of plants explain counterintuitive positive drought effects on an eastern Mediterranean terrestrial ecosystems under climate change. Science 361, 920–923. shrub species. Oikos 119, 1601–1609. Oksanen, J., Blanchet, F.G., Kindt, R., Legendre, P., Minchin, P.R., O’hara, R.B., Simpson, Smith, M.D., Knapp, A.K., Collins, S.L., 2009. A framework for assessing ecosystem G.L., Solymos, P., Stevens, M.H.H., Wagner, H., 2017. Package ‘vegan’. Community dynamics in response to chronic resource alterations induced by global change. ecology package, R package version. Ecology 79, 3–24. Parra, A., Moreno, J.M., 2018. Drought differentially affects the post-fire dynamics of Stevens-Rumann, C.S., Kemp, K.B., Higuera, P.E., Harvey, B.J., Rother, M.T., Donato, D. seeders and resprouters in a Mediterranean shrubland. Sci. Total Environ. 626, C., Morgan, P., Veblen, T.T., 2018. Evidence for declining forest resilience to 1219–1229. wildfires under climate change. Ecol. Lett. 21, 243–252. Pausas, J.G., Llovet, J., Rodrigo, A., Vallejo, R., 2008. Are wildfires a disaster in the Tielbörger, K., Bilton, M.C., Metz, J., Kigel, J., Holzapfel, C., Lebrija-Trejos, E., Mediterranean basin? - A review. Int. J. Wildland Fire 17, 713. Konsens, I., Parag, H.A., Sternberg, M., 2014. Middle-Eastern plant communities Peñuelas, J., Filella, I., Zhang, X., Llorens, L., Ogaya, R., Francisco, L., Pere, C., Marc, E., tolerate 9 years of drought in a multi-site climate manipulation experiment. Nat. Jaume, T., 2004. Complex spatiotemporal phenological shifts as a response to Commun. 5, 5102. rainfall changes. New Phytol. 161, 837–846. Trisos, C.H., Merow, C., Pigot, A.L., 2020. The projected timing of abrupt ecological Peñuelas, J., Sardans, J., Filella, I., Estiarte, M., Llusià, J., Ogaya, R., Carnicer, J., disruption from climate change. Nature 580, 496–501. Bartrons, M., Rivas-Ubach, A., Grau, O., et al., 2018. Assessment of the impacts of van Blerk, J.J., West, A.G., Altwegg, R., Hoffman, M.T., 2021a. Post-fire summer rainfall climate change on Mediterranean terrestrial ecosystems based on data from field differentially affects reseeder and resprouter population recovery in fire-prone experiments and long-term monitored field gradients in Catalonia. Environ. Exp. Bot. shrublands of South Africa. Sci. Total Environ. 788, 147699. 152, 49–59. van Blerk, J.J., West, A.G., Altwegg, R., Hoffman, M.T., 2021b. Does a trade-off between Peterson, B.G., Carl, P., Boudt, K., Bennett, R., Ulrich, J., Zivot, E., 2019. Package growth plasticity and resource conservatism mediate post-fire shrubland responses “PerformanceAnalytics” R Team Cooperation. to rainfall seasonality? New Phytol. 230, 1407–1420. Ploughe, L.W., Jacobs, E.M., Frank, G.S., Greenler, S.M., Smith, M.D., Dukes, J.S., 2019. Wood, S.N., 2017. Generalized Additive Models: an introduction with R, second ed. CRC Community Response to Extreme Drought (CRED): a framework for drought-induced Press. shifts in plant–plant interactions. New Phytol. 222, 52–69. Wright, A.J., Wardle, D.A., Callaway, R., Gaxiola, A., 2017. The overlooked role of Prieto, P., Peñuelas, J., Llusià, J., Asensio, D., Estiarte, M., 2009. Effects of experimental facilitation in biodiversity experiments. Trends Ecol. Evol. 32, 383–390. warming and drought on biomass accumulation in a Mediterranean shrubland. Plant Wu, Z., Dijkstra, P., Koch, G.W., Peñuelas, J., Hungate, B.A., 2011. Responses of Ecol. 205, 179–191. terrestrial ecosystems to temperature and precipitation change: A meta-analysis of Prieto, P., Peñuelas, J., Ogaya, R., Estiarte, M., 2008. Precipitation-dependent flowering experimental manipulation. Glob. Change Biol. 17, 927–942. of Globularia alypum and Erica multiflora in Mediterranean shrubland under Yu, K., Smith, W.K., Trugman, A.T., Condit, R., Hubbell, S.P., Sardans, J., Peng, C., experimental drought and warming, and its inter-annual variability. Ann. Bot. 102, Zhu, K., Peñuelas, J., Cailleret, M., et al., 2019. Pervasive decreases in living 275–285. vegetation carbon turnover time across forest climate zones. PNAS 116, Pugh, T.A.M., Lindeskog, M., Smith, B., Poulter, B., Arneth, A., Haverd, V., Calle, L., 24662–24667. 2019. Role of forest regrowth in global carbon sink dynamics. PNAS 116, Zhang, C., Filella, I., Liu, D., Ogaya, R., Llusià, J., Asensio, D., Peñuelas, J., 2017. 4382–4387. Photochemical Reflectance Index (PRI) for detecting responses of diurnal and Schwalm, C.R., Anderegg, W.R.L., Michalak, A.M., Fisher, J.B., Biondi, F., Koch, G., seasonal photosynthetic activity to experimental drought and warming in a Litvak, M., Ogle, K., Shaw, J.D., Wolf, A., et al., 2017. Global patterns of drought Mediterranean shrubland. Remote Sensing 9. recovery. Nature 548, 202–205. 10