Fire Recurrence and Time Since Last Fire Interact to Determine the Supply of Multiple Ecosystem Services by Mediterranean Forests

Wildfires shape the composition and functioning of Mediterranean ecosystems, but we do not know how these ecosystems respond to both the higher fire recurrence and shorter recovery times expected for future climatic scenarios. We sampled 29 plots with different fire recurrences (from 0 to 4 fires over the past decades) and time since the last fire (up to 35 years; hereafter TSLF) in Southeast Spain, to assess the effect of fire recurrence and TSLF on 25 ecosystem attributes, five related ecosystem services (biodiversity conservation, carbon sequestration, disturbance regulation, food production, and supporting services), plus the synergies and trade-offs between them. High fire recurrence (number of fires) and TSLF interacted to determine ecosystem services but did not affect the synergies and trade-offs between them. Fire recurrence reduced many ecosystem functions and ecosystem multifunctionality. However, this effect dampened, and even became positive, for biodiversity conservation and food production services provided enough (> 20 years) time to recover. The combined effects of fire recurrence and TSLF, however, reduced carbon sequestration and had no overall effects on supporting services. Disturbance regulation, in turn, diminished drastically with the first fire, with no effect of further fires or their interaction with TSLF. Our results show which ecosystem services will suffer more from an increase in fire recurrence, and where restoration and management efforts should focus to maximize the provision of those services more demanded by stakeholders.


INTRODUCTION
Fire affects the composition, structure, and function of vegetation in many fire-prone ecosystems (for example, Mediterranean environments) and increasingly affects ecosystems seldom exposed to fire (for example, tropical forests and temperate deciduous forests) (Mayor and others 2016;Pellegrini and others 2018). Whereas wildfires are a natural disturbance in many forests worldwide, fire recurrence (number of fires experienced over short periods of time) has increased in the last decades (Kovats and others 2014), and it is expected to increase further under future climatic scenarios (Lindner and others 2010;Enright and others 2015). In part of the Mediterranean Basin, for example, fire frequency has doubled over the last century, increasing the burned area by an order of magnitude (Pausas and Ferná ndez-Muñ oz 2012). The effects of these changing fire regimes on vegetation composition, and particularly those of high fire recurrence, are well-known in the Mediterranean region (Santana and others 2010). However, the functional consequences of these compositional changes are less well understood, and even less so the interaction between more recurrent fires and the smaller time windows for ecosystems to recover after these disturbances (that is, time between two subsequent fires; Mayor and others 2016;Karavani and others 2018). Thus, considering TSLF, fire recurrence, and their interaction is important as a forest suffering more fire disturbances is less likely to fully recover (Turner and others 2019). In addition, the degree of recovery is likely to be influenced by how much time the forest had to recover from the last disturbance, and contrasting fire recurrences may alter these post-fire recovery trajectories (for example, Brown and Johnstone 2011;Prichard and others 2017).
In the Mediterranean Basin, many post-fire communities are dominated by non-resprouter seeder species, as a consequence of the long history of land exploitation and subsequent abandonment (Baeza and others 2007). These communities accumulate fine dry biomass and, therefore, can burn again easily (Baeza and others 2006). In these communities, increasing fire recurrence can be associated with higher species richness caused by more open canopies and available niches (Delitti and others 2005). However, this positive effect can dampen as a result of immaturity risk; that is, if fires are frequent enough to prevent plant species to reach reproductive stages and therefore halt their regeneration or lead to local species extinction (Keeley and others 1999;Santana and others 2014). Shorter fire intervals, in addition, may interact with lower post-fire regeneration capacity induced by climate change with drier conditions, leading to a syndrome termed ''interval squeeze'' (Enright and others 2015). This syndrome results in a simplification of vegetation structure that affects other ecosystem functions and services, for example, a reduction in habitat provision for birds and mammals. The lack of trees hinders seed dispersal by birds or mammals, therefore delaying colonization of woody resprouting species (Pausas and others 2006), which concomitantly could prevent vegetation's recovery after further fires. In addition, the lower productivity often observed in frequently burned ecosystems could reduce their ability to sequester carbon (Ló pez-Poma and Bautista 2014; Cruz-Alonso and others 2019). In contrast, the dominance of leguminous shrubs often seen after fires could enhance soil fertility, and more open canopies could favor forage production and understory diversity (Certini 2005;Pausas and Keeley 2019).
In addition to changes in vegetation, fires affect soil nutrients which are crucial for the post-fire recovery of soil biota and vegetation (Duguy and others 2007). Recurrent fires in short periods, through erosion and repeated heating, may progressively produce a long-lasting impoverishment of soil fertility (Knicker 2007). Conversely, fire recurrence can enhance soil carbon and nutrients' concentration by promoting the establishment of more productive plant species, as well as by ash inputs into soils (Boerner and others 2009). Critically, however, studies that document changes in soils affected by more than one fire offer limited insight into long-term changes (Pereira and others 2018). The latter hinders our ability to predict recovery time for these properties and to predict ecosystem service supply under different disturbance scenarios.
Due to the multiple and contrasting effects expected for fire recurrence and TSLF on multiple ecosystem attributes, we need to study not only responses of plant composition but also of the multiple ecosystem functions and services they supply. Furthermore, the different ecosystem functions, and the services supported by them, are not independent of each other, as there may exist synergies (for example, forage productivity and carbon sequestration) or trade-offs (for example, cultural services and timber or food production are difficult to maximize simultaneously; Felipe-Lucia and others 2018). By affecting different ecosystem functions in contrasting ways, fire recurrence and/ Fire Recurrence and Time Since Last Fire Interact to… or TSLF may drive synergies and trade-offs between the ecosystem services to which they relate. The effect of management in these synergies and tradeoffs has received some attention (Raudsepp-Hearne and others 2010; Felipe-Lucia and others 2018), but the effect of natural disturbances has not (Turner and others 2013). In fact, the few studies assessing the relationship between forest fires and multiple ecosystem services focus on few empirical measurements (Pausas and Keeley 2019) and very few of them focus on Mediterranean forests (but see Lucas-Borja and others 2020), which are among the ecosystems most affected by wildfires.
We assessed the effects of fire recurrence (5 levels: 0, 1, 2, 3, and 4 fires), TSLF (from 1 to up to 35 years after the last fire), and their combined effect, on ecosystem multifunctionality (the ability of ecosystems to deliver multiple ecosystem functions and services at high levels simultaneously). We also assessed the effects of fire recurrence and TSLF on five ecosystem services (biodiversity conservation, carbon sequestration, disturbance regulation, food production, and supporting services), and on 25 above-and belowground ecosystem attributes individually. We asked the following questions: (i) how does fire recurrence affect the supply of ecosystem services and overall multifunctionality? (ii) does TSLF interact with fire recurrence to determine ecosystem services? and (iii) do fire recurrence and TSLF affect the synergies and trade-offs between ecosystem services?

Study Sites
The study was conducted in 29 plots (30 9 20 m) located in seven sites in Valencia Region, southeast Spain (Table S1). All plots were old crop fields abandoned in the last half of the twentieth century. These fields have been colonized by Pinus halepensis, with shrubs dominating those plots where fire recurrence was faster than the age of trees to reach maturity. These pine forests and their respective fire-affected stages are representative of natural and semi-natural systems that dominate extensive areas over the Mediterranean Basin (Sheffer 2012 Three of the seven sites were burned experimentally (five recurrences, one fire/10 years, and an unburned control; see Santana and others, 2010). The remaining four sites were affected by natural fires (> 500 Has in surface) with the same range in fire recurrence (number of fires over the last 35 years), but contrasting TSLF. We mixed observational and experimental studies to (i) minimize the correlation between fire recurrence and TSLF ( Figure S1), (ii) double our sample size, and (iii) enlarge the geographic representation of our study. Because the experimental fires were performed in late spring or early summer and under meteorological conditions close to summer wildfires, they reached levels of severity similar to large natural fires (70-90% of biomass lost to combustion, soil surface temperatures up to 710°C, fire intensity between 1600 and 2500 kW m -1 ; Baeza and others 2002; Santana and others 2011). Indeed, sensitivity analyses revealed that the mix between observational and experimental studies did not importantly affect our results and conclusions as we found no significant interactions between the ''study approach'' and our main predictors (fire recurrence and TSLF) for most of our response variables (excepting topsoil organic C, pH, and electrical conductivity; Tables S2 and S3). Mean time-interval between the first and the last fire in our study was 13 ± 5 years for areas burned twice, and 26 ± 5 years for areas burned 3 or 4 times (Table S1). These fire intervals are representative of our study area (Valencia region), where the mean fire return interval for areas burned twice is 14 years. For areas burned 3 or more times, the mean fire interval between the first and the last fire is 29 years, with a gap between fires of approximately 10 years for areas burned 3 times, and of 5-6 years for those with more than 3 fires. Therefore, our experimental design is depicting the fire regime shift toward shorter fire return intervals occurring in the last half-century in Mediterranean areas (Pausas and Ferná ndez-Muñ oz 2012).

Field Sampling and Laboratory Analysis
In each plot, 25 above-and belowground attributes related to their biodiversity conservation potential, their ability to capture carbon on woody biomass and the soil, their potential resistance to additional disturbances (either further fires or soil erosion), their potential to produce food (forage for livestock or wild animals, potential production of honey), or their ability to capture, store and recycle nutrients were collected (see Table 1, Table S4 and text below for details).
In each plot, five 10-m-long transects were set to assess plant species cover, diversity (species richness, Shannon diversity index) using the pointintersect method (Greig-Smith 1983). Measurements were taken every 20 cm and across five heights (0-5 cm, 5-15 cm, 15-25 cm, 25-50 cm, 50-150 cm) (250 points/plot). Potential beekeeping was estimated as the abundance of honeyproducing plants, according to published literature (Sanchís and others 1992;Mateu 2016). The density of woody resprouting species recruiting at each plot was assessed in a subplot of 10 9 10 m that was randomly displayed within each plot.
Shrub and herbaceous biomass were measured in five 1 9 1 m quadrats per plot (one quadrat/transect), and all vegetation and litter were harvested. In the laboratory, the clipped material was separated into dead and live biomass. Samples were then oven-dried and weighed. The biomass of large woody species was measured by allometric equations using the basal stem diameter (Baeza and Santana 2015) of all tree individuals found in a 10 9 10 m square (the same subplot used to measure the density of resprouter species). The total woody biomass was calculated by summing shrub and tree biomass. Dead shrub biomass and the vertical shrub cover were used to calculate the dead fuel bulk density (fuel weight / volume; hereafter DFBD) as a surrogate for flammability (that is dead shrub biomass (g/m 2 ) / shrub height (m); (Santana and others 2011).
To analyze the effect of fire recurrence on soil attributes, five soil samples (in the same five 1 9 1 m quadrats as the biomass samples) within each study plot were collected. All soils were sampled during early summer, to avoid seasonal variation in the variables measured. The soil samples were air-dried and sieved (< 2 mm) before the analyses. Soil pH and electrical conductivity (hereafter EC) were determined in a mixture of  (Nelson and Sommers 1996). Available phosphorus was obtained following Olsen and others (1954) method. The enzymatic activities bglucosidase and acid phosphatase related to the cycling of carbon and phosphorus were analyzed following the procedure described by Tabatabai (1994). Each soil attribute, except soil enzymatic activities, was measured into two soil subsamples at two depths (0-5 and 5-15 cm). As an additional assessment of the functional status of each study plot, the landscape functional analysis (LFA; Tongway and Hindley 2004) was used to measure the stability, infiltration, and nutrient cycling indices (Appendix S1). After removing those highly correlated with others, and therefore providing redundant information (Table 1; Figure S2), the different vegetation and soil variables were organized, by calculating their averages after being standardized between 0-1, into five categories of ecosystem services: (i) carbon sequestration: live woody biomass and soil organic carbon, as the main sources of stable carbon storage in terrestrial ecosystems (Heimann and Reichstein 2008), (ii) biodiversity conservation: plant species richness (related also to the diversity of insects and belowground biota; Scherber and others 2010), and habitat complexity (a composite index related to habitat provision for birds, lizards and mammals, using the cover of different vegetation strata and litter, according to Val and others (2017), (iii) food production: beekeeping potential and herbaceous biomass (as an indicator of forage availability for both livestock and game), both with a high economic importance in the study area, (iv) disturbance regulation: related to the ability of the ecosystem to resist soil erosion (stability index), likelihood of getting burned again (DFBD, multiplied by -1 to reflect resistance to further fires), and the ability to recover after an additional fire (resprouter abundance), and (v) supporting: nutrient cycling index, b-glucosidase, and acid phosphatase. We calculated a multifunctionality index as the mean of the five ecosystem services (see Allan and others 2015 for a related approach). Alternatively, instead of the averaging approach, we also calculated multifunctionality using the multiple threshold function defined as the number of functions reaching some percentage (40, 60 and 80% in our case) of the maximum observed (Byrnes and others 2014). Since the results from the latter were very similar to those presented here ( Figure S3), we do not discuss them further.
Finally, we analyzed the trade-offs and synergies between ecosystem services applying the method used by Felipe-Lucia and others (2018). To do so, first, we calculated the Pearson correlations between the standardized ecosystem services (raw synergies and trade-offs). Second, we removed the effect of the combined effect of fire recurrence and TSLF by taking the residuals from the linear mixed model (Fire recurrence and TSLF as predictors) and then calculating the Pearson correlations between ecosystem services again (intrinsic synergies and trade-offs, regardless of fire regime). Third, we removed the effect of fire recurrence (using the residuals with it, but not TSLF, as a predictor; recalculating the correlations again; synergies and trade-offs related to fire recurrence). Fourth, we removed the effect of TSLF (using residuals with it, but not fire recurrence, and recalculating the correlations between the ecosystem services again: synergies and trade-offs related to TSLF). Finally, we evaluated differences between the four correlation matrices (raw data, residuals after fire recurrence, residuals after TSLF, and residuals after both fire recurrence and TSLF), to assess differences in the synergies (positive correlations) and tradeoffs (negative correlations) between the ecosystem services evaluated (see Appendix S2 for more details).

Statistical Analyses
To assess the effect of fire recurrence and TSLF on each ecosystem function and service, we performed a linear mixed-effects models analysis using the R function lme from the package nlme (Pinheiro and others 2020) and the function r.squaredGLMM from the MuMIn to calculate the R-square for the mixed effect models (Burnham and Anderson 2002). Adequacy of model fits to our data, and lack of bias was assessed using the ''EnvStats'' package. To explore which relationship fits better our data, we compared linear vs quadratic relationships for ecosystem services and multifunctionality and found that the linear relationship was the most parsimonious model fitting our data (Table S5). We used fire recurrence and TSLF, and their interaction, as fixed factors, and study site as a random factor (random intercept). Before this, we removed the influence of the different climatic conditions characterizing each study site by fitting a model with each of our response variables and the mean temperature and rainfall of each site (as predictors) and taking the residuals for further analyses. Before any analysis, all variables were standardized between 0 and 1 to obtain comparable coefficients.
We could not include the forests with 0 fires in these analyses, as they do not have an associated TSLF. Thus, we only used them as a standard to compare the rest of the sites, and we did so by comparing the level of multifunctionality and the five services between 0 fires and the rest of fire recurrences with a mixed model using ''fire recurrence'' as the sole fixed predictor, and site as a random factor. We also re-calculated each value of ecosystem services and multifunctionality obtained for the different fire regimes as a proportion of the average of the unburned forests, to aid interpretation of these unitless metrics. All statistical analyses were performed with R version 3.5.2 (R Core Team 2017).

Ecosystem Multifunctionality, Services, and Functions
Increasing fire recurrence reduced ecosystem multifunctionality (t = -2.28; P = 0.04, Figures 1, 2 and S5; Table S6 and S8). Multifunctionality was 20% lower (regarding unburned forests) in forests burned 1 to 3 times, and 70% lower in forests burned four times. Fire recurrence and TSLF interacted to determine three out of the five services measured: biodiversity conservation, food production, and carbon sequestration (Figure 1; Table S6). In two cases (biodiversity conservation and food production), the supply declined with fire recurrence, although this negative effect dampened, and even shifted to positive if TSLF were large enough (20 years or higher). Carbon sequestration responded positively to TSLF (t = 2.41; P = 0.03), but declined with fire recurrence, particularly after the second fire (matching tree loss in most cases). The response of carbon sequestration was driven mainly by changes in woody biomass (t = -2.86; P = 0.01, Figure S4; Table S7). The disturbance regulation and supporting services did not respond to changing fire regimes. However, when compared with the unburned plots, these two ecosystem services declined after the first fire by ca. 40% and 30%, respectively ( Figure 2; Table S8). In addition, the effect of fire recurrence on dead fuel bulk density, and important attribute related to further fire risk (and therefore associated with disturbance regulation services), was negative (t = -2.97; P = 0.01, Figure S4; Table S7).
Overall, the interaction between fire recurrence and TSLF, rather than their individual effects, affected most ecosystem services and aboveground functions ( Figures. 1 and S4; Tables S5 and S7), while no significant effects were observed on belowground variables ( Figure S6; Table S9).

Trade-offs and Synergies between Ecosystem Services
Correlations between pairs of ecosystem services ranged between -0.11 to 0.68 (Figure 3a), and only two negative correlations were observed: those between disturbance regulation and carbon sequestration ( -0.11) and between carbon sequestration and biodiversity regulation ( -0.03). When removing the effect of TSLF (Figure 3c), we observed a (non-significant) decline in most of the synergies observed between ecosystem services, with some of them even shifting to trade-offs (between carbon sequestration and food production; Tables S10 and S11). We did not find significant changes in the relationships between ecosystem services when removing the effect of fire recurrence, TSLF, or their combined effect. The latter suggests that the synergies and the trade-offs observed are intrinsic to the supply of these services, and not driven by the compositional, structural, or functional changes we observed in response to fire recurrence and TSLF. Although these results could be also driven by our relatively small sample size, perhaps not big enough to statistically determine this influence, they did not change when including the six extra unburned sites in the analyses (data not shown). In any case, we found strong and consistent synergies between three of the five ecosystem services studied in all cases: biodiversity conservation, food production, and disturbance regulation ( Figure 3).

DISCUSSION
Our study shows that changes in fire regimes can substantially alter the supply of important ecosystem services and functions from Mediterranean forests. Ecosystem multifunctionality, that is, the number of ecosystem services supplied simultaneously at high levels, declined by a factor of three when fire recurrence increased from one to four fires (and by 70% regarding unburned forests), regardless of the time allowed to recover. If wildfires become more frequent, as expected (Turco and others 2014), this would critically impair the ability of Mediterranean forests to store carbon, halt biodiversity loss, or produce food. These negative effects on Mediterranean ecosystem services may be particularly severe in the context of climate change, where squeezed fire intervals are expected Fire Recurrence and Time Since Last Fire Interact to… (Enright and others 2015) and therefore would not allow these ecosystems enough time to recover. Ecosystem services and functions responded in different ways to the effects of fire recurrence and TSLF (Figures 1; S4 and S6; Tables S6; S7 and S9). Some, such as disturbance regulation (risk to further fires, ability to resprout after a disturbance, resistance to soil erosion), showed strong declines after the first fire to never fully recover, whereas others, such as soil attributes, showed no responses overall to both fire recurrence and TSLF (Figures 2  and S6; Tables S8 and S9). Our results suggest, therefore, that ecosystems experiencing wildfires more frequently or for the first time, such as some boreal or tropical forests (Aragã o and others 2018; Whitman and others 2019), are unlikely to fully recover during the next decades (Brown and Johnstone 2011).
Importantly, long TSLF may buffer the negative effects of fire recurrence and even enhance the supply (positive interaction between fire recurrence and TSLF) of other ecosystem services such as biodiversity conservation and food production.
Thus, in order to enhance the supply of these services, management and surveillance should aim to minimize fire hazards in those ecosystems already burned several times, to allow them enough time to recover. The latter could be achieved either by reducing fuel load or re-introducing species with low flammability (for example, resprouter species) in frequently burned areas (Santana and others 2018). It must be noted, however, that these scenarios (high fire recurrence followed by long TSLF) are extrapolations outside of our data range, and have not been observed in any case for fires larger than 500 Has in the entire study region since records exist ( Figure S1). Furthermore, if recurrent fires affect a given site before trees have had the opportunity to reach the adult stage and refill their aerial seed bank, this may produce drastic declines in tree biomass, ecosystem productivity, and resilience to further disturbances (Tapias and others 2001;Mayor and others 2016;Pellegrini and others 2018;Turner and others 2019), strongly reducing its ability to sequester carbon (negative fire recurrence 9 TSLF interactions; Figure 1; Table S6) Similarly, in northern boreal and Yellowstones lodgepole pine forests, carbon stocks dramatically declined following short-interval fires (Brown and Johnstone 2011;Turner and others 2019). Such vegetation changes may be catastrophic in terms of carbon sequestration and climate change mitigation (Santana and others 2016; Turner and others 2019), and our results extend those previous findings to the ability of the ecosystem to resist soil erosion or further fires.
Despite the strong declines in ecosystem multifunctionality, and also on many of the aboveground attributes measured, fire recurrence and TSLF did not strongly affect soil attributes (Table S9). This result is partially at odds with findings in forests with more developed soil organic layers (Pellegrini and others 2020) and could be related to the lack of measurements related to N availability, one of the most sensitive soil nutrients to fires, in our study. However, results in other Mediterranean forests have also reported similarly weak effects of fires on soils, and N availability in our study region seems not to respond strongly to changes in fire regime (Moghli and others unpublished data; Gué non and others 2013). Conversely, the lack of effects could also be partly due to our mixture of observational and experimental data. Indeed, including the study approach (experimental vs natural fires) in our analyses induced changes in our results regarding some of the belowground attributes, with stronger negative effects on soils found on the natural versus experimental fires (Table S3). These few changes (mainly in the topsoil) are probably related to fire duration (probably higher in natural fires than in the experimental ones) may be the main factor that caused the greatest belowground differences in natural versus experimental fires (Certini 2005). Regardless of the former, the generally weak influence of changing fire regimes on soil attributes is consistent with previous findings in cork oak woodlands, where soil carbon was similar between different fire regimes (1-2 vs. 4 fires) after 17 years of recovery (Gué non and others 2013). Similarly, Ferran and others (2005) did not find significant differences between soil organic carbon at 2.5-10 cm when comparing fire recurrences.
Our study shows that fire recurrence may reduce both the number of species (particularly when fire recurrence > 3) and the variation in functional traits (by removing trees), unless ecosystems are allowed to recover for more than 20 years. Provided enough time, however, canopy openings and the mix between late and early successional species may have a positive effect, as shown by the positive effect of the interactions between fire recurrence and TSLF. This interaction suggests that repeated fires reduce initially the richness of species as a consequence of seed bank depletion (Goudelis and others 2007; Santana and others 2014), but long Values higher or lower than 100% indicate that burned forest supply more or less of that particular ecosystem service than the unburned forest reference, respectively. To ease the visualization, the values of biodiversity conservation within plots burned twice (225.3%) and three times (176.5%) are not shown in the graph, and values belonging to each fire recurrence are slightly scattered across the X axis but separated by vertical dashed lines.
Fire Recurrence and Time Since Last Fire Interact to… TSLF may allow again the entrance of opportunistic and heliophyte species in the opened gaps. In contrast, in plots with low fire recurrence, a long TSLF reduces the diversity due to secondary succession (that is, after the tree canopy closes, species richness drops reaching its lowest values 35 years after a fire, or in unburned forests vs those burned 1 or 2 times; Figure 2). Furthermore, many studies found that the spatial and temporal heterogeneity of vegetation structure may provide a greater range of habitat for animals (Fuhlendorf and others 2006;Val and others 2017). We found that the combined effect of fire recurrence and TSLF reduces this habitat complexity ( Figure S4). This result is a consequence of tree and litter removal, and it is consistent with previously observed declines in heterogeneity after recurrent fires (for example, in tallgrass prairies, Collins 1992). All these changes in diversity may have large consequences for overall ecosystem multifunctionality provided the oft-reported positive biodiversity-functioning relationships in woodlands (Jucker and others 2014; Liang and others 2016; but see van der Plas 2019). In-deed, we found a strong correlation between biodiversity conservation and overall ecosystem multifunctionality (r = 0.81, P < 0.001), which is explained by the positive relationship between biodiversity and forest functioning (for example, Jucker and others 2014), the positive effect of heterogeneous forest structures on ecosystem multifunctionality (Felipe-Lucia and others 2018), and the positive linkages between biodiversity conservation and other ecosystem services (food production, disturbance regulation, and supporting services Figure 3). In addition, our results highlighted the generally positive effect of fire recurrence on the food production attributes we studied. Herbaceous biomass increases with the combined effect of fire recurrence and TSLF, whereas the woodland shifts to shrubland may increase the beekeeping potential (that is, strong correlation (r >|0.7| observed between shrub biomass and beekeeping potential). The positive responses of herbaceous layers to fires are associated with nutrient flushes into the soil early after the fire and the presence of open canopies. These are well- known effects, as fire has been historically used to enhance forage production for livestock and game in many places (Pausas and Keeley 2019). Given their contrasting responses, future studies should address the effect of fire recurrence and TSLF on other food production attributes, such as game, mushroom, wild berries, or edible plants (see Gamfeldt and others 2013;Gassibe and others 2014;Felipe-Lucia and others 2018).
We did not find a significant effect of fire recurrence or TSLF on disturbance regulation or supporting services. This can be explained by the contrasting responses of the ecosystem attributes feeding these services, or by the weak responses of most soil attributes (our supporting services) to fire recurrence and TSLF, as discussed above. Resprouters' abundance declined with fire recurrence, whereas dead fuel biomass responded to the combined effect of fire recurrence and TSLF (decreased fire risk with recurrent fires and long TSLF), and those soil surface indices related to stability did not respond at all. The latter could be related to the long-timescale covered by our study, as strong rates of soil erosion are commonly observed during the first months after the fire, but those decline rapidly with vegetation recovery (Cerda 1998).
Our study showed that fire recurrence, TSLF, and their combined effect do not induce shifts in the trade-offs and synergies commonly found between ecosystem services. This result suggests that the relationships among ecosystem services provided by Mediterranean forests are consistent and not driven by the strong structural and compositional changes in vegetation we observed under differing fire regimes (Figure 3, Tables S10 and S11). Our results are partially at odds with those found in temperate forests, as the synergies and trade-offs of the services they provide were strongly affected by management-related structural changes (Felipe-Lucia and others 2018). However, our results agree partly with these previous findings, as biodiversity conservation (related to habitat heterogeneity) was positively correlated with food production and disturbance regulation services across all our sites. The lack of response of ecosystem services' tradeoffs and synergies to fire recurrence and TSLF could be explained by two main causes: (i) the prevalence of soil measurements in our study (which respond more weakly to these structural changes) or (ii) the influence of landscape features that were not considered here. The latter could alter the provision of some services linked to connectivity between different ecosystem types (for example, biodiversity conservation) or spatial location (for example, erosion risk is the highest upslope), and therefore modulate the synergies and trade-offs in between these services. Although non-significant, our results suggest that fire recurrence may reduce the synergies between ecosystem services or even change them to trade-offs. This result warrants further research to assess whether TSLF and the associated structural and compositional changes are an important factor enhancing ecosystem synergies, as suggested both by Felipe-Lucia and others (2018) in temperate forests and by Lucas-Borja and Delgado-Baquerizo (2019) in Mediterranean ecosystems.

CONCLUSIONS
Our study shows the main changes in five important ecosystem services supplied by Mediterranean forests, the synergies, and trade-offs between them, and multifunctionality overall, as affected by high fire recurrence and TSLF. We illustrated the importance of considering different timescales when studying the effects of fire recurrence, as they interact strongly to determine ecosystemic responses. Our study also highlights the need to consider the diversity of functions and services that natural ecosystems supply, as they respond differently to fire recurrence (and recover at a different pace). This gives important information for land managers and policymakers as to which ecosystems services can be expected under different fire regimes, and which require strong fire prevention or active restoration efforts to keep them at high levels (for example, disturbance regulation). Our findings suggest that, when managing at the landscape level, it is very important to maintain a vegetation mosaic to maximize the supply of multiple ecosystem services simultaneously. Thus, sustaining pine woodland may increase carbon sequestration and the colonization of resprouters may enhance disturbance regulation. However, the use of fire for promoting shrub and open patches could be useful for enhancing other ecosystem services. For example, repeated prescribed burning combined with a subsequent long-term fire exclusion can enhance services such as food production and diversity, as well as reduce fuel loads at the landscape level. We also found strong synergies between ecosystem services related to biodiversity conservation, food production, and disturbance regulation which would pave the way to manage multifunctional Mediterranean landscapes that satisfy the needs of stakeholders with little ecological costs.
Fire Recurrence and Time Since Last Fire Interact to… AC KNOWLEDGEMENT S A. Garrapiso, F. Romero, M. Drakopoulou, L. Stout, and V. Ferná ndez helped with soil analyses and during fieldwork. J.A. Alloza helped with plot selection and the calculation of the fire return interval. We thank the anonymous reviewers for their constructive comments that improved the paper. This research was funded by the Spanish Society of Terrestrial Ecology and the generosity of A.Traveset (AEET; project REFRESCOME, awarded to S.S.), the HYDROMED project funded by the Spanish Ministry of Science and Innovation (Subprojects BLUEWATER PID2019-111332RB-C21 and INERTIA PID2019-111332RB-C22), and the FIRE-SCENARIO project within the program for emerging research groups funded by the Generalitat Valenciana (GV-2020-160). A.M. is supported by the scholarship of Generalitat Valenciana-European Social Fund (ACIF-2018-194). S.S. and E.P. were supported by the Spanish Government under a Ramó n y Cajal contract (RYC-2016-20604).

FUNDING
Open Access funding provided thanks to the CRUE-CSIC agreement with Springer Nature.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/li censes/by/4.0/.