**Projecting Future Impacts of Global Change Including Fires on Soil Erosion to Anticipate Better Land Management in the Forests of NW Portugal**

**Amandine Valérie Pastor 1,2,\*, Joao Pedro Nunes 2, Rossano Ciampalini 1, Myke Koopmans 3, Jantiene Baartman 3, Frédéric Huard 4, Tomas Calheiros 2, Yves Le-Bissonnais 1, Jan Jacob Keizer <sup>5</sup> and Damien Raclot <sup>1</sup>**


Received: 1 October 2019; Accepted: 4 December 2019; Published: 11 December 2019

**Abstract:** Wildfire is known to create the pre-conditions leading to accelerated soil erosion. Unfortunately, its occurrence is expected to increase with climate change. The objective of this study was to assess the impacts of fire on runoff and soil erosion in a context of global change, and to evaluate the effectiveness of mulching as a post-fire erosion mitigation measure. For this, the long-term soil erosion model LandSoil was calibrated for a Mediterranean catchment in north-central Portugal that burnt in 2011. LandSoil was then applied for a 20-year period to quantify the separate and combined hydrological and erosion impacts of fire frequency and of post-fire mulching using four plausible site-specific land use and management scenarios (S1. business-as-usual, S2. market-oriented, S3. environmental protection and S4. sustainable trade-off) and an intermediate climate change scenario Representative Concentration Pathway (RCP) 4.5 by 2050. The obtained results showed that: (i) fire had a reduced impact on runoff generation in the studied catchment (<5%) but a marked impact on sediment yield (*SY*) by about 30%; (ii) eucalypt intensification combined with climate change and fires can increase *SY* by threefold and (iii) post-fire mulching, combined with riparian vegetation maintenance/restoration and reduced tillage at the landscape level, was highly effective to mitigate soil erosion under global change and associated, increased fire frequency (up to 50% reduction). This study shows how field monitoring data can be combined with numerical erosion modeling to segregate the prominent processes occurring in post forest fire conditions and find the best management pathways to meet international goals on achieving land degradation neutrality (LDN).

**Keywords:** sediment yield; runoff; fire frequency; erosion control techniques; mulching; global change

#### **1. Introduction**

Erosion is a serious environmental problem worldwide including on-site effects (e.g., depletion in soil organic carbon stocks and decline in agronomic yields) and off-site effects (e.g., non-point source pollutions and reservoir siltation) [1]. Although it is a natural phenomenon, it can be exacerbated by climate change; fire and grazing [2] and is often accelerated by human activities [3]. The problem is

especially disruptive in the Mediterranean area where erosion rates measured at the outlet of river basins of different sizes are high [4,5], especially because soils are usually thinner than in northern Europe [6] but also because croplands of the Mediterranean have a predominance of fallow land and lack of adoption of conservation practices [7]. The Mediterranean basin is very prone to erosion in all its forms because of the climate, topography, soil characteristics and a very long history of human presence and intense cultivation [8,9]. In addition, the occurrence of human-related forest fires affecting thousands of hectares each year is also a significant problem in both the northern and southern areas of the Mediterranean basin because wildfires have frequently been reported to produce strong and sometimes extreme hydrological and erosion responses in recently burnt areas, especially during the first few post-fire years [10,11].

Mediterranean-type ecosystems are particularly prone to fires, as the cool wet season is propitious for vegetation growth and fuel production, while the summer dry season promotes the drying of fuel and propitious conditions for fire occurrence and spread [12]. In the north-western Iberian Peninsula, the more humid climate has favored the commercial plantation of fast-growing but fire-prone species such as eucalypts and maritime pines, further exacerbating the wildfire problem [13]. Besides damaging vegetation and leading to loss of life and property, fires have a large number of indirect impacts, which may be more important in the long term, such as impacts on health through smoke inhalation; disruptions of social functions such as road and air traffic or business closures; environmental impacts in burned areas such as enhanced soil erosion, flooding and water contamination in and downstream of burned areas and long-term disruption of social-ecological dynamics [14]. In the Mediterranean, the enhancement of soil erosion in burned areas is particularly disruptive because, once the forest soils are directly exposed to the water action, they can be rapidly degraded by post-fire erosion and associated carbon and nutrient losses [15,16]. The export of fine sediments, associated with ashes and nutrients, can also contaminate streams and impact aquatic ecosystems and human water resources [17–19].

Fires are strongly dependent on climate [12,20]. Therefore, climate change projections of warmer and, in many fire-prone regions, drier climates indicate an increase of fire frequency and extend in the Mediterranean [21] and other regions [22,23]. However, fires are also dependent on vegetation, both in terms of potential vegetation growth and fuel production [12] and of changes to vegetation distribution, caused by climate or socioeconomic changes [24,25]. The resulting changes in wildfire regimes may lead to changes in their impacts, and requires an assessment of both the magnitude of these changes and the effectiveness of existing measures to mitigate them. The impacts of wildfire on soil fertility losses can often be effectively mitigated through soil conservation measures, especially when applied as emergency post-fire interventions. These measures include the application of organic residues to the soil surface (i.e., mulching, with straw or forest logging residues), log and shrub barriers, infiltration enhancement through scarifying or plowing or re-seeding and ecological restoration [17,26–29]. However, it is typically impractical to treat all burnt areas entirely, thus requiring knowledge of the effectiveness of these measures under different fire and post-fire conditions [29–34].

The impact of fire frequency on post-fire hydrological and erosion response has not yet been fully established for the forests of NW Iberia (Portugal and Spain) [35], and approaches to prioritize post-fire recovery are not widely applied in practice, except in Galicia [34]. More knowledge is also needed about the potential combined impacts of fires and global change on runoff and soil erosion, given the potential to exacerbate the already damaging impacts of fires on land degradation. The objective of this study was to assess the potential of soil conservation measures, including post-fire mulching, to mitigate the impacts of global change in a Mediterranean forested catchment of NW Portugal.

To this end, the long-term soil erosion model LandSoil was first calibrated using rainfall-runoff-erosion measurements on a catchment that burnt in 2011, and then applied over a 20-year period to quantify the individual and combined impacts on runoff and sediment yield (*SY*) of fire frequency (from no fire to two fires every 20 years), post-fire mulching and global change scenarios derived from four plausible land use and management scenarios under an intermediate climate change scenario.

#### **2. Material and Methods**

#### *2.1. Study Site Macieira*

The study was carried out in the Macieira de Alcôba experimental watershed, of 96 ha, located at an altitude of 470 m asl in the Caramulo mountain range of north-central Portugal (Figure 1). The climate can be classified as humid Mediterranean (Csb in the Koppen-Geiger classification) with annual precipitations ranging from 818 to 1294 mm (Pousadas rainfall station). About 70% of the precipitation falls during autumn and winter. The average air temperature in winter is 8 ◦C, while that in summer is 17 ◦C. Soils are shallow, mostly Luvisols and Cambisols, overlaying moderately impermeable bedrock of schist and granite [15,19]. The Southern part of the watershed contains the village surrounded by terraced agricultural fields with a rotation of maize in summer and pasture in winter, irrigated year-round with the "águas da lima" system [19]. Until the 1930s, the watershed was mainly occupied by agricultural fields on terraces and by shrublands on the steepest slopes; from the 1930s to 1970s, there was a large-scale afforestation of shrublands and part of the agricultural lands with Maritime Pine (*Pinus pinaster*); from the 1970s to present, afforestation with eucalypt (*Eucalyptus globulus*) occurred [36]. Around the village and agricultural area, the forest is present on slopes steeper than 5%. With a reliable map of burnt areas starting in the mid-1970s, major fires have affected the region about once per decade [37].

**Figure 1.** Location of Macieira de Alcoba, Portugal.

#### *2.2. Data Collection*

The Macieira watershed was monitored between 01 November 2010 and 30 September 2014, with continuous measurements of rainfall and other climatic variables, and streamflow and water turbidity at the outlet (Table 1). Observation showed that sediment exportation at the outlet is mainly due to suspended material because bedload is negligible in this site [28]. A fire burned c. 10% of the watershed on 11th August 2011, and about half of the burnt area was clear-cut and ploughed in February and March 2012; therefore, the available data set included two periods: before the fire (2010–2011) and after fire (2012–2014). Outlet data collection was complemented by a survey of linear erosion features composed of rills and ephemeral gullies after major storms, which the resulting volumes of eroded soil converted to mass using bulk density measurements. Furthermore, a characterization of soil texture, depth, hydrological properties (water retention curve, saturated hydraulic conductivity and others), roughness and shear strength, was carried out through a survey at 25 points within the catchment. A more detailed description of data collection can be found in Nunes et al. (2018) [38].


**Table 1.** Geomorpological and Hydrological Research Unit (HRU) characterization.

#### *2.3. Model Setup and Calibration Procedure*

The LandSoil model (LandSoil: landscape design for soil conservation under land use and climate change) is a spatial raster-based model for simulating water and tillage erosion and landscape topography evolution over time spans on the orders of 10–100 years [39,40]. After simulating each event, LandSoil recalculates the elevation raster. Soil surface properties control water infiltration, runoff and sediment concentration for each grid-cell and rainfall event. As the result of water infiltration balance, runoff is routed over the catchment using a D8 algorithm [41] integrating the upstream/upslope cells contribution and the runoff generated within each cell. Linear landscape elements as ditches and tillage rows are also used to impose preferential directions to the flow [42].

LandSoil can simulate both rill and inter-rill erosion processes. The first is based on factors such as runoff friction and cohesion, topography (slope) and accumulated flow in an empirical relationship with the observed rill sections [43]. The latter, responsible for the remobilization of the soil particles detached by splash erosion, is controlled by the sediment concentration responding to rainfall intensity and soil surface properties such as vegetation cover, soil roughness and soil crusting as declined in expert rules issue of observations at parcel scale [44,45]. In LandSoil, sediment concentration is limited by threshold functions combining local topography, including profile curvature (concavity > 0.055 m−<sup>l</sup> ) and slope gradient (<0.02 mm−), and land use/vegetation cover (>60%). Sediment concentration limits the range from 2.5 to 10 gl−. Tillage erosion is simulated following the formulation of Govers et al. (1994) [46], both contour and downslope were modeled in different fields with tillage transport coefficients spanning from 111 to 139 kgm−<sup>g</sup> [47].

The implementation of LandSoil on the studied catchment was based on a 2 m resolution grid based on a digital elevation model (DEM) derived from interpolating 1:10,000 contour lines using the ridge line method [48]. The set-up for LandSoil modeling mainly consists on defining the monthly calendars for soil roughness and crusting, vegetation cover, and tillage operation per land use based on a field sampling campaign in 2010 to determine soil texture, chemistry and origin (Table 1, Tables S4–S9). A digital elevation model (DEM) derived from interpolating 1:10,000 contour lines using the ridge line method [48], at 25 m of spatial resolution has been used to represent the catchment topography. In LandSoil, climate forcing is provided through a time series of precipitation events expressed through the effective volume and duration of precipitation, previous precipitation over 48 h and maximum precipitation over five minutes. The calibration of the model was made using events that occurred during the four years monitoring period in a 2-steps procedure: first for runoff, and then for erosion. Calibration for runoff consisted in manually adjusting the soil infiltration parameter. Calibration for erosion consisted in manually adjusting the sediment concentration in runoff and the rill sections. Parameters optimization was based on runoff and sediment exportation at the catchment outlet, as well as field measurements of hydrological soil properties and rill distribution sizes. In LandSoil, soil and land use parameters evolution are provided at a monthly base although simulations are made at the event rainfall base. This choice is consistent with the objective of simulating runoff and erosion

over time pans of more than 10 years as it ensures a good representativity of erosion at the catchment scale [49].

All the events generating effective runoff and corresponding to a duration above 3 h, volume above 40 mm and maximum 5-min intensity above 2 mm/h were selected and used in the calibration process. This selection includes both pre- and post-fire events generating together about 90% of the total measured soil erosion. Additionally, in order to inspect a priori the influence of the rainfall factors on runoff and erosion, a regression analysis was performed to relate runoff and *SY* with rainfall volume, antecedent rainfall, maximum intensity and duration of the event.

#### *2.4. Climate and Fire-Related Scenarios*

In this study, climate data were downloaded from the platform Medcordex.eu for both historical (1986–2005) and future (2041–2060) periods. Only the intermediate Representative Concentration Pathway (RCP) 4.5 climate emission scenario [50] was used because differences between RCP scenarios by 2050 were very small. The Regional Climate Model (RCM) ALADIN was chosen in the Mediterranean Agricultural Soils Conservation under global Change (MASCC) project; it was considered as a valid climate model in the comparison of RCM by Fantini et al. (2018) [51]. A statistical downscaling based on a quantile mapping approach was applied on the daily time rainfall series for bias correction, and a temporal downscaling based on an analog approach was then applied to transform the daily rainfall time series into five minutes' time series. Thus, the 5-min time series have a real temporal structure as a direct result of the observations.

Burnt area scenarios were calculated from climate scenarios according to Sousa et al. (2015) [21]. To this end, the annual burnt area in the region surrounding Macieira (the Águeda river basin) between 1990–2016 was calculated from Portuguese historical databases [13], and a multiple linear regression was calculated with monthly rainfall anomalies, extreme fire weather (days with daily severity rating (DSR) above the 90th percentile [52]) and fire history in the preceding years. Meteorological data was taken from the ECMWF Re-Analysis [53]. While the agreement with observations was not very high (r2 = 0.50), it was sufficient to provide a "best guess" of how much and when the burned area would change for the climate scenarios. Burnt area scenarios under climate change were calculated using rainfall and DSR projections, the latter with a correction for bias using an empirical quantile mapping approach [54]. The estimated doubling of burnt area, reflecting an increase in fire frequency rather than the size of the burnt areas, was more conservative than the three- to four-fold increase obtained with the same methodology for northwestern Iberia [21]. The results for the Águeda were subsequently downscaled for Macieira watershed, by assuming a doubling of fire frequency, from 1× to 2× every 20 years. Fire location was that of the 2011 fire, while the burnt area was the same in all scenarios, except for S2 with 5% larger burn area due to the eucalypt stand in the southern part of the catchment (Figure 2).

#### *2.5. Land Use and Management Scenarios*

Four narrative storylines and land use scenarios were used in this study (Figures 3 and 4, Supplementary S1). For this, six local researchers working in different institutions of the Central Portugal NUTS II region (where Macieira is located) were interviewed and completed a survey on writing narrative scenario storylines for four contrasted socio-economics trajectories, initially, designated in the MASCC project: S1 (business-as-usual), S2 (market-oriented), S3 (environmental protection) and S4 (sustainable) in the context of global change in the Caramulo region (Supplementary S1, Tables S1–S3). The survey encompassed a table with future land use distribution, a table with allocation constraints per land use, and a ranking exercise on post-fire soil management rated by local experts. The main elements of the narratives of the four expert-based scenarios are shown in Figure 3.

**Figure 2.** Location of the burned area in orange in our simulations corresponding to the area that burned in 2011.

**Figure 3.** The four narrative land use scenarios elaborated under the MASCC project and used in this study.

**Figure 4.** Land use scenarios.

The current situation (S0) in Macieira is described as that of traditional agriculture (including corn and grassland) that has gradually been abandoned for afforestation, especially with *Eucalyptus* plantations that now occupy both abandoned agricultural areas and large areas of former pine forest and shrubland (Figure 4). Orchards, although still marginal, are increasing in all our scenarios as being a new profitable production system (berries). In all the scenarios, pine plantation decreases due to low profitability and due to reported increased pests and diseases, while oak is re-introduced in the environmental protection scenario (S3) and the sustainable scenario (S4) for biodiversity conservation and increase in soil water holding capacity, with incentives from local associations and/or government. In the S3 scenario, the maintenance and/or restoration of riparian vegetation is foreseen to capture sediments along the stream (with vegetated strips), in combination with adoption of conservation agriculture with lower tillage intensity. In S4, the same biodiversity and soil protection measures were implemented as in S3 but to only half the extent. In S3 and S4, soil tillage was also adjusted to follow contour lines.

As post-fire management, we selected application of mulch because it was suggested by local-experts as the most cost-effective post-fire conservation measure in this region, which is consistent with the findings of studies in the region [28,29,32,55] as well as in other parts of Iberian Peninsula [56,57].

#### *2.6. Simulation Design*

Model simulations concerned the period from 1986 to 2005 as well as a future period of 20 years between 2041 and 2060. They included the simulation of the exhaustive time series of major rainfall events, i.e., rainfall events with a mean lager than 40 mm. First, we assessed the impact of fire frequency and mulch application on runoff and soil erosion (Table 2) and then we assessed the combined impacts of land use and climate change scenarios, following the matrix in (Table 3).

**Table 2.** Impact of fire frequency and mulching on runoff and sediment yield (*SY*) on the baseline scenario S0.

> **No Climate Change (noCC)** No fire 1 fire 1 fire + Mulch 2 fires 2 fires + Mulch


**Table 3.** Impact of global change on runoff and sediment yield (*SY*).

#### *2.7. Data Analysis*

The impact of fire frequency, mulch application and changes in climate and land use and management on runoff and erosion were analyzed by comparing the LANDSOIL simulations outputs with R studio software Version 1.0.136. In this paper, three levels of fire frequency were tested (i.e., 0, 1 and 2 fires every 20 years); two levels of mulch application (with or without application); two levels of climate change (with or without) and five scenarios of land use and management (S0 to S4). The simulation of the combinations of factors described in Tables 3 and 4 permitted to evaluate the distinct and combined impacts of these factors for the 2041–2060 period compared to the 1986–2005 period. A period of 3 post-fire years was also used to analyze the direct impact of fire and related contribution over 20 years. When there was no fire, we calculated the runoff and erosion produced during the period when fire would have had occurred to be able to compare with the period with fire. When two fires occurred, we calculated the runoff and *SY* that occurred two times: for three years after each fire.

**Table 4.** Impact of fire frequency and mulch application on runoff and sediment yield (*SY*) over 20 years (1986–2005). Control represents «S0 + 1 Fire» in 20 y, ii represents «no fire», iii represents «1 Fire including post-fire mulch application», iv represents «2 Fires» and v represents «2 Fires including post-fire mulch application». Relative values are in italics.


\* Fire events represents the 3 years after the fire occurred.

#### **3. Results**

#### *3.1. Variable Explaining Runo*ff *and Sediment Yield in Macieira de Alcoba*

Multiple linear regression analyses of event-wise runoff volumes (in m3) and *SY*s (in Mg) gave satisfactory fits. This was especially true for runoff with an R2 of 0.85 as opposed to 0.58 for *SY*. The best-fitting equations were:

$$Runoff = 581.12 \times Pt - 97.28 \times PI \\ 30 - 165.74 \times Pd + 213.34 \times AP \\ 48 - 20354.82, \tag{1}$$

$$SY = 49.46 \times Pt - 18.55 \times PI \\ 30 - 14.64 \times Pd + 131.97 \times AP \\ 48 - 2964.59,$$

where: *Pt* is rainfall volume (in mm), *PI30* is the maximum 30-min rainfall intensity (in mm.h<sup>−</sup>1), *Pd* is the duration of the rainfall event and *AP48* is the amount of rainfall over the 48 h preceding the event. The regressions results showed that rainfall and *SY* were explained best by *Pt* (*p* < 0.001), then by *AP48* (*p* < 0.01) and less significantly *Pd* and *PI30* (see Supplementary S1).

#### *3.2. Calibration of Runo*ff *and Sediment Yields Using the LandSoil Model (2010–2014)*

The calibration results were obtained from creating two sets of infiltration rates according to the rainfall events (pre-fire and 3 years of post-fire disturbance) and by creating two sets of erosion by calibrating sediment concentration for diffuse erosion and rill dimension according to pre-fire conditions and fire conditions (three years following fire event) based on the observed data. The calibration of 30 events resulted in a satisfactory simulated runoff with R<sup>2</sup> of 0.92 and R2 of 0.85 for simulated *SY* (Figure 5). We observed a better performance of the simulated runoff and sediment yield in the period post-fire than pre-fire however the period before fire had much lower values than after fire. The averages and SD were similar in both periods and variables. The NSE has the total value of 0.91 for runoff and 0.86 for sediment yield (see Supplementary S9).

**Figure 5.** Results of calibration of LandSoil model for runoff (**left plot**) and sediment yield (**right plot**) with observed events between 2010 and 2013. SD stands for sediment yield.

#### *3.3. Impact of Climate Change on Major Precipitation Events*

Under the future climate conditions, the total number of major rainfall events over the 20-year period was 5% higher than under present climate conditions. Furthermore, the cumulative rainfall of these events was 7% higher under the future than present climate conditions (Figure 6 and Supplementary S1). At the same time, rainfall duration and maximum intensity, as well as antecedent rainfall, were between 15% and 30% higher under future climate conditions than under the current situation.

**Figure 6.** Boxplots of key characteristics of major rainfall events under present (no climate change—noCC, 1986–2005) and future (climate change—CC, 2041–2060) climate conditions: rainfall total (**a**,**b**), maximum intensity over 5 min (**c**,**d**), duration (**e**,**f**) and antecedent 48 h rainfall (**g**,**h**).

#### *3.4. Impacts of Fire Frequency and Post-Fire Mulching on Erosion under Present Climate Conditions*

Compared to the control scenario (one fire over the 20 years period 1986–2005), the scenario without fire reduced runoff by 3% over 20 years and by 23% during the 3 years of post-fire disturbance while the scenario with two fires increased runoff by 5% over 20 years and by 153% during the 2 × 3 years of post-fire disturbance (Table 4). In turn, post-fire mulching reduced runoff by 1% and 2% in the case of one and two fires over 20 years, respectively and up to 7% during the post-fire period only. Sediment yield was more strongly affected by fire frequency and post-fire erosion mitigation than runoff. The scenario without fire decreased *SY* by 33% over the 20 year period (and by 82% during the post-fire period only) compared to the control scenario, while the one with two fires increased it by 28% over the 20 years and by 86% (during the post-fire period only). The scenarios involving post-fire mulch application decreased *SY* by 16% and 28% in the case of one and two fires, respectively (and by 35% and 40% during the post-fire period only). While with one fire, the years of fire disturbance represents 15% of the runoff generation, when doubling fire frequency, the two periods of fire disturbance generation about 37% of runoff. Concerning *SY*, under the control situation, the fire disturbance period represents 43% of the *SY* but when fire frequency doubles, 69% of the sediments are produced during both fire periods (six years out of 20). In the scenario without fire, we showed that during the 3 years when fire would have occurred; only 12% and 15% of runoff and *SY* were respectively produced. Finally, the contribution of interrill erosion to total erosion accounted for 1%–8% of the total erosion in all the simulations.

#### *3.5. Runo*ff *and Sediment Yield Under Future Climate Conditions*

Without changes in fire frequency, cumulative runoff under future climate conditions was 9%–21% higher than under present climate conditions (Table 5, Figures 7 and 8). With a doubling of fire frequency, the cumulative runoff was 13%–25% higher under the future than present climate conditions. The impacts of changes in climate conditions, as well as fire frequency, were most pronounced in the case of the business-as-usual (S1) and market-oriented scenarios (S2). Land use had reduced impacts on runoff under present climate conditions (with differences of up to 3%, in the case of S2) but marked impacts under future climate conditions, with runoff differences amounting to 11%–13% for the S3 and S4 scenarios and to 22%–23% for the S1 and S2 scenarios.

**Table 5.** Impact of fire frequency under land use and management and climate change on cumulative runoff and sediment yield (*SY*) over 20 years. Runoff is in mm and sediment yield in Mg ha<sup>−</sup>1. The impact is quantified as a percentage of change relative to the control. noCC represents no climate change within the 1986–2005 period and CC represents climate change within the 2041–2060 period.


**Figure 7.** Cumulative sediment yield (*SY*) in Mg ha−<sup>1</sup> vs. rainfall in mm between 1986 and 2005 with one fire occurrence in 1995.

**Figure 8.** Cumulative sediment yield (*SY*) in Mg ha−<sup>1</sup> vs. rainfall in mm under climate change between 2041 and 2060 with two fire occurrences in 1941 and 1951.

Concerning sediment yield (Table 5), we observed that the impact of climate change increased *SY* from 11% (S4 scenario) up to 45% in the S1 scenario under one fire and the impact of increasing fire frequency could increase *SY* up to 29%. Concerning the impact of land use, we observed high contrasts between scenarios. Without climate change, *SY* increased from 32% in S1 up to 90% in S2 while we observed a decrease in *SY* in S3 and S4, respectively by 52% and 21%. The combined impacts of climate change and land use change show an increase of *SY* by 60% and 126% in S1 and S2 respectively and a decrease by 22% and 47% in S4 and S3 respectively. Land use change and management have a strong impact on *SY* with nearly a doubling of *SY* in S2 and a reduction by half in S3 scenario. Figure 9 shows the net cumulative erosion (soil loss) and deposition (soil gain) between 1986 and 2005 for the S2 scenario. Three fields, including two burned fields, show the highest cumulative soil loss between 7 and 9 Mg ha−1, while the areas in the north of the catchment show net deposition probably coming from the road located on the NW of the catchment. In Figure 10, the impact of increasing riparian vegetation in S3 shows the impact of sediment deposition around the stream. Overall, we show that the highest negative impact would be to change the landscape to mainly *Eucalyptus* plantations with no conservative land management approaches (S2); while with sustainable management (S3) *SY* can be reduced by half.

**Figure 9.** Soil loss between 1986 and 2005 in S2 after two fires (Kg). Negative values represent soil loss and positive values represent gain.

**Figure 10.** Soil amount difference between Scenario S3 and S0 without a fire at year 2005 (Kg). Negative values represent loss and positive values represent gain.

#### **4. Discussion**

#### *4.1. Application of LandSoil to Macieira*

In our study, we found that the volume of rainfall during and 48 h before the fire event was the most important climatic driver of the catchment-scale soil erosion. From these, rainfall volume had the strongest correlation with erosion rates and this can be due to the humid Mediterranean climate in Macieira de Alcoba, with usually bigger and longer events than in drier Mediterranean regions, and concurs with the results of similar other studies in NW Iberia [58]. The calibration of the infiltration rate, sediment concentration and rill dimension led to a good performance of LandSoil in representing current runoff and *SY*. The only land use that had to be calibrated according to literature was the re-introduction of oak forest, which leads in time to a higher soil water infiltration than Eucalypt [59]). Performance of LandSoil might be improved with an extended field measurement dataset but this would have required more field monitoring and probably another fire occurrence with a different extent and location in the catchment.

#### *4.2. Overall LandSoil Results*

The model application to Macieira indicates that the impact of climate change is relatively low, increasing runoff generation by 9% and sediment yield by 19%. Doubling fire frequency has low impacts on runoff in the long-term, increasing by 5%, but larger impacts on *SY,* which increases by 28%. However, the model also shows the effectiveness of mulch application after fires on *SY* with a reduction of 16%–28% with one and two fires respectively (and up to 40% reduction during the post-fire disturbance period only). The application of mulch had a lesser impact on runoff than on *SY* (up 2% reduction) over 20 years and up to 7% during the post-fire disturbance period. Finally, the model indicates that land use change had the lowest impact on runoff, with changes between −1% and 3%, but the largest impacts on *SY*, changing between −23% and 92%. These results show that land use allocation and management was the most important factor impacting *SY* with a potential exacerbation of fire and climate change impacts in the case of S1 and S2 scenarios and strong potential of mitigation in S3 and S4 scenarios. In these latter scenarios, a large part of sediments was trapped in riparian fields, therefore preventing sediments to reach the stream and outlet (Figure 10).

This decreased *SY* more than soil erosion; while this would not mitigate land degradation in burned areas, it would limit stream contamination by ash transport after a fire event, which can have significant impacts on water resources [19]. Compared to the modeling study by Nunes et al. (2018b), our results showed less impact of fire on runoff and *SY*; while in Nunes et al. (2018a) [60] a 'no fire situation' decreases *SY* by 2/3, in our case *SY* was decreased by 1/3. This could be due to the duration of the calculated cumulative *SY*. In Nunes et al. (2018a, b), this is done for a 10-yr period as opposed to for a 20 years period in this study. In fact, when our results were only compared during the three years of post-fire disturbance, the efficiency of mulch application decreased runoff by 7% and *SY* by 40%. Compared to Prats et al. (2019, 2012) where mulch application reduced runoff by 15%–25% and *SY* by more than five times, our results were more conservative probably because our burnt field area only represented 10% of the watershed area; while Prats et al. (2012, 2019) calculated *SY* at the field area scale and they burned 100% of the area while in our case 10% of the catchment burned. Moreover, our simulations ran over 20 years vs. 5–10 years in other studies [32,61].

Our study shows the low impact of mulch application on runoff generation. This could be due to high soil water repellency in the eucalypt burnt forest. Soil water repellency also exists in *Eucalyptus* forests and in a pre-fire situation due to preponderant ligneous root materials compared to oak and pine trees [59,62]. However, soil water repellency has a limited impact on runoff generation, in the early autumn [61]; and this is only important for patches since at the hillslope scale there are always non-repellent zones where water can re-infiltrate [15]. The low impact of mulching on runoff reduction is also probably due to the low total amount of runoff generation in both burnt and unburnt eucalypt plantations, at least at the hillslope scale [15,37]. Finally, the risk of increasing floods remains

relatively low compared to predicted extreme rainfall events [60] due to the low impact of fires on runoff production.

Concerning sediment reduction caused by mulching, Shakesby et al. (1996) [63], showed that *Eucalyptus* logging litter reduced soil erosion by 95%. The average reduction was five times in Prats et al. (2019, 2012); it should be noted that, in the latter study, maximum erosion rates were lower than in our study site. As previously said, the study site of Prats et al. (2012) also presents smaller temporal and spatial scales. In our study, mulch could reduce *SY* by up to 40%, which was quite high considering that fire only impacted roughly 10% of the watershed. We also suggest that these values vary with the type and quantity of mulch that is applied, the soil type, the percentage of soil cover and the slope of the site (Prats et al., 2014, 2012). Finally, Figures 5 and 6 show that erosion was mainly a raindrop-driven process with the detachment of sediments in eucalypt forest and some tilled field. The highest levels of soil loss occurred in the burnt field area and in one field with a slope above 10◦ in the S2 scenario (about 7–8 Mg per parcel) and these three areas contribute to 37% of total soil loss in scenario S2 (Figures 7 and 8). The combined impact of fire on *SY* at the outlet and the spatial representation of soil loss coming from burned areas show that erosion was mainly due to soil detachment rather than transported by water to the outlet via sedimentation on the riparian field close to the stream and in the terraces.

#### *4.3. Distinct and Combined Impacts of Global Change and Fire Frequency*

In light of future climate change, it is likely that fire frequency could increase in a context of climate change and that doubling fire frequency might be conservative in our study when compared with burned area projections for this region [21]. The land use scenarios might also impact fire frequency and extent; for example, the increase of *Eucalyptus* plantation and intensification in S1 and S2 might exacerbate the risk of fire occurrence and propagation, while in S3 and S4 the adoption of crop diversification and protective land use and management should lead to lower risks of fires and erosion. However, we did not take these differences into account because we wanted to have the same number of fires in our land use scenarios so the specific impact of fire would be comparable. In any case, our results show that fire frequency had a larger impact on soil erosion (up to 33% change) than on runoff generation (up to 5% change). In addition, our results also show 43% of the erosion was produced during 3 years after the fire event representing 15% of the duration of the simulation (Table 4). This concurred with multiple studies, which showed that the indirect impacts of climate change, such as associated changes to land use and fire occurrence can have larger impacts than those of changes to climatic variables alone [9,64–66]. In the fire-prone forests of NW Iberia, efforts to limit the impact of climate change on land degradation should, therefore, focus on fire prevention and post-fire impact mitigation.

#### *4.4. Status of Post Fire Measures Implementation*

Our study shows the potential of using soil management measures to mitigate post-fire erosion in NW Portugal with a mix of numerical modeling and field experiments. The biophysical effectiveness of erosion control measures is well understood in NW Portugal [67,68]. Studies show the high performance of mulching and erosion control barriers, and the lesser effectiveness of erosion control measures such as post-fire seedling, which are more difficult to manage due to the time lag of seed germination and the availability of native seeds [69]. However, there is currently only one study calculating the costs of post-fire control erosion measures [68]. The benefits produced by controlling overland runoff and *SY* with control erosion measures on hydrological services should be biophysically and economically quantified so that implementation of erosion control measures can be communicated at the science-policy management level [70]. For instance, the risk of water contamination must be studied and costs of prevention (with erosion control measures) compared to the costs of water treatment must be evaluated [18,60,71,72]. The same applies to the hydrological service of flood regulation and enhanced water soil content and fertility with post-fire erosion control measures via

increased soil organic matter, soil aggregate stability and nutrient availability [73]. For this work to proceed, more studies are needed on the socio-economic aspects of fire prevention and post-fire contamination control.

#### **5. Conclusions**

This study evaluated the potential impacts of global change and fire frequency on water runoff and soil erosion. The novelty of this study lies with the assessment of global change with field and landscape conservation management techniques such as application of mulch in post-fire conditions, implementation of vegetated strips (or restoration of riparian vegetation) and decreased soil tillage on water runoff and soil erosion. Our results show that the LandSoil model performed well for runoff and *SY* compared to the observed data of 2010–2014 and that rainfall volumes (>40 mm) during and 48 h before fire occurrence are the main climate drivers of soil erosion before rainfall duration and maximum intensity. With LandSoil we could test contrasted land use scenarios including *Eucalyptus* intensification to show their potential risk of implementation under global change (up to +98% in *SY* in the market-oriented (S2) scenario) while a conservation scenario including environmental protection and afforestation of local species can decrease soil erosion by half. The application of mulching was an effective tool against erosion (up to 28% with two fires in 20 years). Finally, preventing fire is also a potential erosion mitigation tool with a decrease of 33% of soil erosion. Valuating hydrological services of post-fire erosion measures would clarify the feasibility of implementation of these latter by communicating the results within a science-policy interface in a context of the latest Intergovernmental Panel on Climate Change (IPCC) recommendations on achieving climate change mitigation and land degradation neutrality (LDN).

**Supplementary Materials:** The following are available online at http://www.mdpi.com/2073-4441/11/12/2617/s1, Supplementary S1: Narrative Storyline of Land Use Scenario of MASCC Project, Supplementary S2: Land Use and Climate Scenarios. Table S1: Land use allocation constrains for Macieira, Table S2: Land use allocation per scenario S0 (baseline 2000), S1, S2, S3 and S4, Table S3. Soil conservation practices ranking by six local experts during MASCC project, Table S4: Statistical summary of rainfall events selected on the study site during the period 1986–2005 (noCC) and 2041–2060 (CC), Supplementary S3: Results from Statistical Analyses of Observed Rainfall Events and Corresponding Runoff and Sediment Yield (R version 3.4.2 (28 September 2017)—"Short Summary"), Supplementary S4: Table S4: Calibration of Monthly Calendars of Soil Roughness, Soil Crusting and Soil Cover per Land Use, Supplementary S5: Soil Tillage Calibration, Supplementary S6: LandSoil Turbidity Calibration, Supplementary S7: Calibration of Rill Erosion, Supplementary S8: Calibration of Infiltration Rates for Two Periods (Pre-Fire and during the Fire Disturbance Period), Supplementary S9: Results from the Calibration of Runoff and Sediment Yield in Pre and Post Fire Periods. Figure S1: Difference in elevation change between scenario S3 and S0 without fire occurrence. Figure S2: Elevation change in scenario S2 between 2000 and 2050.

**Author Contributions:** Conceptualization: A.V.P., J.P.N., R.C., Y.L.-B., J.J.K. and D.R.; Formal analysis: A.V.P., J.P.N., M.K. and F.H.; Funding acquisition: T.C. and D.R.; Investigation: J.P.N., Y.L.-B., J.J.K. and D.R.; Methodology, A.V.P., J.B., F.H., T.C., Y.L.-B. and D.R.; Software: R.C. and M.K.; Supervision: J.P.N., J.B. and Y.L.-B.; Validation: J.P.N. and D.R.; Visualization: A.V.P.; Writing—original draft: A.V.P., J.P.N., D.R.; Writing—review & editing, J.P.N., R.C., M.K., J.B., F.H., T.C., J.J.K. and D.R.

**Funding:** This work benefits from the financial support of the MASCC project through ARIMNET2, an ERA-NET funded by the European Union s Seventh Framework Program for research, technological development and demonstration under grant agreement no. 618127. APV was funded by the MASCC project described above. This work was further funded by the "Fundação para a Ciência e a Tecnologia", with personal grants attributed to J.P. Nunes (IF/00586/2015) and J.J. Keizer (IF/01465/2015), and funding for data collection in Macieira within project ERLAND (FCOMP-01-0124-FEDER-008534). J. Baartman visited the site with funding from an STSM mission within the Connecteur COST Action ES1306. M. Koopmans was partly funded by ERASMUS+.

**Acknowledgments:** We would like to thank the two external reviewers and the editorial board for their valuable comments before the acceptance of the manuscript.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

1. Lal, R. Accelerated soil erosion as a source of atmospheric CO2. *Soil Tillage Res.* **2019**, *188*, 35–40. [CrossRef]


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

## *Article* **Multiple Temporal Scales Assessment in the Hydrological Response of Small Mediterranean-Climate Catchments**

**Josep Fortesa 1,2,\*, Jérôme Latron 3, Julián García-Comendador 1,2, Miquel Tomàs-Burguera 4, Jaume Company 1,2, Aleix Calsamiglia 1,2 and Joan Estrany 1,2**


Received: 13 November 2019; Accepted: 14 January 2020; Published: 19 January 2020

**Abstract:** Mediterranean-climate catchments are characterized by significant spatial and temporal hydrological variability caused by the interaction of natural as well human-induced abiotic and biotic factors. This study investigates the non-linearity of rainfall-runoff relationship at multiple temporal scales in representative small Mediterranean-climate catchments (i.e., <10 km2) to achieve a better understanding of their hydrological response. The rainfall-runoff relationship was evaluated in 43 catchments at annual and event—203 events in 12 of these 43 catchments—scales. A linear rainfall-runoff relationship was observed at an annual scale, with a higher scatter in pervious (R2: 0.47) than impervious catchments (R2: 0.82). Larger scattering was observed at the event scale, although pervious lithology and agricultural land use promoted significant rainfall-runoff linear relations in winter and spring. These relationships were particularly analysed during five hydrological years in the Es Fangar catchment (3.35 km2; Mallorca, Spain) as a temporal downscaling to assess the intra-annual variability, elucidating whether antecedent wetness conditions played a significant role in runoff generation. The assessment of rainfall-runoff relationships under contrasted lithology, land use and seasonality is a useful approach to improve the hydrological modelling of global change scenarios in small catchments where the linearity and non-linearity of the hydrological response—at multiple temporal scales—can inherently co-exist in Mediterranean-climate catchments.

**Keywords:** rainfall-runoff; multiple temporal scales; non-linearity; small catchments; Mediterranean

#### **1. Introduction**

The complexity of Mediterranean fluvial systems is caused by the multiple temporal and spatial heterogeneities in the relationships between the natural and human-induced abiotic and biotic variables, especially in the Mediterranean Sea Region [1,2]. The Mediterranean climate lies between 32◦ and 40◦ N and S of the Equator and is characterized by a wet and mild winter, a warm and dry summer and a high inter- and intra-annual variability in rainfall patterns. Mediterranean climate regions comprise the Mediterranean Sea Region, the coast of California, Central Chile, the Cape region of South Africa and the southwestern and southern parts of Australia [3]. Under these climatic conditions, catchments are mostly characterized by a high diversity in hydrological regimes [4,5] promoting significant temporal and spatial differences in the hydrological response [6–8]. The seasonality of the Mediterranean climate plays a key role in the runoff generation processes, increasing the non-linearity of the rainfall-runoff relationship at the event scale [9–11]. In winter and early spring, saturation processes are dominant, due to water reserves triggering the runoff generation [8,12]. During late spring, summer and early autumn, those same authors observed how runoff was generated under Hortonian conditions due to high rainfall intensities. Different runoff mechanisms can co-exist within a catchment [13], although flood events under antecedent saturation wetness conditions enable a larger hydrological response [12,14–16].

The spatial variability of runoff generation is also elucidated in the catchment hydrological response as a combination of rainfall-distribution [17] and runoff-contribution areas [18]. Thus, the spatial distribution of landscape elements in relation to each other is fundamental in influencing transfer flow pathways [19]. The main factors governing connectivity are associated to changes in topography [20], soil properties [21], soil type [22] and in vegetation cover [23]. In addition, river connectivity occurs along a three spatial dimension (i.e., longitudinal, lateral and vertical) which can change markedly over time, especially in ephemeral rivers [24]. Consequently, flow pathway activation depends on rainfall amount and intensity, as well as soil moisture antecedent conditions [25–27]. Lithology's effects on runoff response in Mediterranean Sea Region fluvial systems are conditioned by the presence of karst features, as the proportion of carbonate rocks is significantly higher than in other landscapes [28]. Therefore, its characteristics related to hydrology, such as high infiltration rates, deep percolation and spring sources, must be taken into account [29–31]. In limestone areas, Hortonian and saturation runoff can both be generated and infiltrated downslope [32,33], whilst in badland areas the runoff generation is characterized by a lower soil infiltration capacity than in karst areas [34]. In areas with high clay subsoil content [35] and over granite bedrock [36], runoff is generated as a combination of a lack of deep percolation and a subsurface runoff over an impervious subsoil layer, even artificially promoted [15].

Land uses also alter the hydrological response, depicting a runoff reduction when agriculture uses are replaced by forests [37]. Afforestation processes due to land abandonment can promote a 40% reduction in the annual water yield in Mediterranean Sea Region fluvial systems [38], whilst forest logging may increase the annual runoff coefficient by up to 16% [39]. However, afforested catchments generate the largest flows and peak discharges at the event scale when compared with forest catchments [40]. Most of these catchments are affected by soil and water conservation structures historically built to reduce overland flow and prevent erosion [41]. The abandonment and degradation of these structures may promote and increase runoff and sediment yield [42], with runoff coefficients being between 20% and 40% in abandoned terraces [43].

Hydrologists develop perceptual models of the catchments they study, which consist of an appreciation of the dominant processes controlling the hydrological response, using field measurements and observations [44]. Numerical models used to simulate catchment behaviour often fail build on this knowledge [45], also considering that the Mediterranean-climate catchments show very heterogeneous responses over time and space, resulting in limitations in hydrological modelling and large uncertainties in predictions [46]. To reduce this spatio-temporal scale variability, small experimental and representative catchments are useful to observe the hydrological response under different or specific land use, lithology and human effect characteristics [47]. The aim of this paper is to investigate the rainfall-runoff relationship at different temporal scales in representative small Mediterranean-climate catchments (i.e., <10 km2), evaluating the role of lithology and land use. At the annual scale, the runoff response was assessed at 43 catchments under a pervious or impervious lithology. At the event scale, the rainfall-runoff response of 203 events was investigated to examine the effects of seasonality, lithology and land use. In addition, the inter- and intra-annual variability of the rainfall-runoff and the temporal downscaling (i.e., annual to event scale) was studied in the Es Fangar Creek catchment (3.35 km2; Mallorca, Spain) during five hydrological years. The rainfall-runoff

relationship assessment under contrasted lithology, land use and seasonality may provide new insights into the hydrological modelling of small Mediterranean-climate catchments. These catchments are highly demanding in terms of data and event-scale modelling because classical hydrological models are not well suited to the Mediterranean area, as many hydrological processes remain poorly represented [46].

#### **2. Materials and Methods**

#### *2.1. Study Areas*

#### 2.1.1. Small Mediterranean-Climate Catchments

A total of 43 small catchments (i.e., <10 km2) from 22 published studies on the Mediterranean climate regions (Figure 1a) were selected to analyse hydrological response. The geographical distribution of the catchments was grouped into the main climate regions, as follows: (1) Western coast of USA, (2) Western Mediterranean Sea Region (from Spain to Italy), (3) Eastern Mediterranean Sea Region (Israel) and (4) South Africa. The area of the catchments ranged from 0.05 to 9.61 km2, the median value being 1.03 km2 and the standard deviation 2.6 km2. The mean annual rainfall ranged from 367 to 1794 mm y−1, with a median value of 833 mm y−<sup>1</sup> <sup>±</sup> 334 mm y−1. The mean annual temperature ranged from 6.6 to 17.2 ◦C with a median value of 13.9 ◦C ± 3 ◦C. When temperature information was not available, it was obtained according to Fick and Hijmans [48]. The predominant lithology was pervious in 12 catchments, and was impervious in the other 31. Within the 43 catchments, the studies of 12 of them also contained information related to the main land uses, which was used in this paper for assessing their hydrological response at the event scale. The main land uses were agriculture (3 catchments), agroforestry (3), forestry (1) and shrub (5).

#### 2.1.2. Es Fangar Creek

A temporal downscaling assessment of the inter- and intra-annual variability of the rainfall-runoff relationship was carried out in the Es Fangar Creek catchment, a headwater tributary of the Sant Miquel River catchment (151 km2) located in the north-eastern part of Mallorca Island (Figure 1b) and is representative of Mediterranean mid-mountainous catchments. The lithology is mainly composed of marl and marl-limestone formations from the medium–upper Jurassic and Cretaceous period in the valley bottoms. In the upper parts of the catchment, massive calcareous and dolomite materials from the lower Jurassic period and dolomite and marl formations from the Triassic period (Rhaetian) are dominant. The Es Fangar catchment has an area of 3.4 km2, with altitudes ranging from 72 m.a.s.l. to 404 m.a.s.l. (Figure 1c). The mean slope of the catchment is 26% and the length of the main channel is 3.1 km (average slope of 22%). The drainage network is natural in the headwater parts. In the bottom valley, flow lamination is applied, with transverse walls and also the straightening and diverting of the main stream, with the banks fixed with dry-stone walls for flood control and erosion prevention. In addition, subsurface tile drains are also installed to facilitate drainage due to the impervious materials which would impede agricultural activity during wet periods. As a result, 16% of the surface catchment is occupied by soil and water conservation structures. Since 1950, important socio-economic changes have caused a gradual abandonment of farmland in marginal areas, leading to afforestation. The land uses in 1956 were rainfed herbaceous crops (54%), forest (31%) and scrubland (15%). Nowadays, the main land uses (Figure 1c) are forest (63%), rainfed herbaceous crops (32%) and scrubland (5%). In addition, 54% of terraced land is currently covered by forests (Figure 1c), demonstrating the consolidation of the forest transition. The climate of the area is classified on the Emberger scale [49] as Mediterranean temperate sub-humid. The mean annual rainfall (1965–2016, Biniatró AEMET station) is 927 mm y−<sup>1</sup> with a variation coefficient of 23%, and the mean annual temperature is 15.7 ◦C. A rainfall amount of 180 mm in 24 h is estimated to have a recurrence period of 25 years [50]. The Es Fangar streamflow regime can be classified as intermittent flashy (49% zero days

flow), with an annual variability from intermittent (35% zero days flow) to harsh intermittent (62% zero days flow).

**Figure 1.** (**a**) Map of the small Mediterranean-climate catchments selected to assess the rainfall-runoff relationship at the annual and event scale. (**b**) Map of Mallorca Island, showing the location of the Sant Miquel River and Es Fangar Creek catchments. (**c**) Map of the Es Fangar Creek catchment, showing the different land-uses, the stream network and the gauging station. (**d**) Map of the Sant Miquel River catchment with the location of rainfall stations used in this study.

#### *2.2. Hydrological Response of Small Mediterranean-Climate Catchments*

Bivariate statistical regressions were used to establish the correlations at the annual and event scales between rainfall and runoff in order to assess the hydrological response of small Mediterranean-climate catchments (i.e., <10 km2; Figure 1c). These equations were used descriptively as least squares adjusted for the rainfall-runoff assessment. At the annual scale, data from the 43 representative catchments (Table A1 (a); hereinafter A1a) were collected to observe the influence of lithology on this response; i.e., the catchments were classified as pervious or impervious by using the information regarding the catchments' characteristics (e.g., soil type, soil texture or lithology materials) extracted from research papers. At the event scale, 203 events from 12 representative catchments were classified according to (a) seasonal occurrence (autumn, winter, spring or summer), (b) pervious or impervious lithology and (c) main land use (agricultural, agroforestry, forest or shrub) (Table A1 (b)). The main land uses of each catchment were divided into specific land uses when the information was available (Figure 1). Most of these studies (79%) were located in the Mediterranean Sea Region and 47% of them studied catchments < 1 km2.

#### *2.3. Monitoring and Data Acquisition in Es Fangar*

The rainfall data since 2012 were obtained from the B696 Biniatró AEMET station (Figure 1d; 1 km away from the catchment). In October 2014, a rainfall gauge station (Míner Gran) was installed less than 2.5 km away from the catchment. The Míner Gran rainfall gauge is located 1 m above the ground and connected to a HOBO Pendant® G Data Logger-UA-004-64 that records precipitation at 0.2 mm resolution. A linear regression was established (n = 978; R2: 0.88) for daily rainfall (2014–2017) between the Biniatró and Míner Gran stations to reconstruct rainfall data series from 2012 to 2014 for the Míner Gran station. Due to the lack of temperature data available in the studied catchment, the data of neighbouring AEMET weather stations (i.e., less than 8 km away from the catchment) were used to estimate the catchment's temperature by using the block kriging technique. With this information, the monthly evapotranspiration (i.e., ETo) was estimated using the equation from Hargreaves and Samani [51].

The gauging station of the Es Fangar catchment was built in July 2012. Its cross section is formed by a rectangular broad-crested weir for low water stages to better measure low discharges. The water level was continuously (1 min time step) measured using a pressure sensor Campbell CS451 connected to a CR200X datalogger and average readings were kept every 15 min. The flow velocity was measured during baseflow conditions and flood events using an OTT MF Pro electromagnetic water flow meter. These flow velocity measurements (n = 17) were subsequently used to calibrate the stage–discharge relationship.

#### *2.4. Rainfall-Runo*ff *Relationship Assessment in Es Fangar*

This study is based on data from 5 hydrological years (2012–2013 to 2016–2017). For each runoff event (i.e., when the water stage exceeds the low-flows channel; i.e., 0.036 m<sup>3</sup> s<sup>−</sup>1), a simple hydrograph separation between quickflow and baseflow components was performed through a visual technique based on the breakpoints detected on the logarithmic falling limb of the hydrograph [52]. This hydrograph separation method was used only to characterize the response of the catchment to a rainfall, with no aim of deriving any interpretation in terms of runoff processes. Over a 5-year period, hydrograph separation was conducted on 49 events. Based on the hydrograph separations, the baseflow and quickflow contributions for each hydrological year (October to September) were calculated.

Several variables were derived from the hyetograph and hydrograph from each rainfall-runoff event (Table 1). These variables aimed to characterize the pre-event conditions (antecedent precipitation in 24 h, AP1d; baseflow at the start of the event, Q0) as well as the event characteristics (rainfall depth, Ptot; mean and maximum rainfall intensity, IPmean30 and IPmax30; runoff depth and runoff coefficient, R and Rc; and peak-flow discharge, Qmax). The relationships between these variables were assessed

through the Pearson correlation matrix. A more detailed analysis of rainfall-runoff relationships in the Es Fangar catchment was carried out in a second step to investigate the variability of the hydrological response to similar size rainfall events to assess the rainfall-runoff non-linearity. The widest range of events with a similar Ptot (40–50 mm) generating the highest differences in R response was selected. Accordingly, 7 events with a Ptot range from 41.8 to 49.8 mm but with different antecedent conditions or rainfall dynamics were selected and also classified in relation to wet, dry and transition periods, in accordance with Gallart et al. [53].

**Table 1.** Pre-event and event conditions variables to explain the rainfall-runoff relationship in the Es Fangar Creek catchment.


#### **3. Results**

#### *3.1. Hydrological Response of Small Mediterranean-Climate Catchments*

#### 3.1.1. Annual Scale: Lithology Influence

The annual rainfall ranged from 376 to 1794 mm yr−1, whilst R ranged from 3.7 to 1200 mm within these representative catchments under Mediterranean-climate conditions. The relationship between annual rainfall and R (Figure 2) showed a significant positive linear correlation (R2 = 0.68; p < 0.01). However, some scattering was also apparent in the relationship because when catchments with pervious lithology were not included, the regression increased (R2 = 0.82; p < 0.01; see Figure 2). The lowest annual R values (<10 mm) are related to catchments with annual rainfall < 500 mm yr−<sup>1</sup> and a pervious lithology (Figure 2), located in the Eastern Mediterranean Sea (Lower Jordan River). Besides, catchments with pervious lithology with ca. 1500 mm yr−<sup>1</sup> of annual rainfall showed large differences in their R amount, ranging from 70 to 974 mm.

**Figure 2.** Rainfall-runoff for 43 small Mediterranean-climate catchments at the annual scale. Catchments with pervious lithology are marked with a grey halo. The Es Fangar Creek value is illustrated with a black dot.

#### 3.1.2. Event Scale: Seasonality, Lithology and Land Use Influences

To investigate the variability of the hydrological response at the event scale, the rainfall-runoff relationships of 203 events from 12 small Mediterranean-climate catchments were analysed (Table A1 (b) and Figure 3). Within Figure 3, the outliers are marked with an ellipsoid and excluded from the correlations because outliers enlarged the rating curves of these correlations, obtaining significant relationships in all cases.

**Figure 3.** Rainfall-runoff relationship at the event scale classified by (**a**) season, (**b**) lithology and (**c**) land uses at 12 small Mediterranean-climate catchments. Outliers are marked with an ellipsoid.

The seasonal distribution of the events was winter (34%), autumn (32%), spring (22%) and summer (12%)—most of them occurred between November and February (47%). Similar seasonal median rainfall was observed, ranging from 26.2 (winter) to 29.0 mm (autumn). The highest seasonal median of event R occurred in winter (4.2 mm). However, those events occurred during the transition periods depicted a similar median R (autumn = 2.3 mm; spring = 2.2 mm) being the lowest value for summer events (0.6 mm). In the case of rainfall-runoff relationships, the highest correlation was obtained in spring (R<sup>2</sup> = 0.63), followed by winter, autumn and summer (Figure 3a).

Small differences between the median of R events in catchments with pervious lithology (2.2 mm) and catchments without pervious lithology (3.1 mm) were observed. Nevertheless, significant differences in rainfall-runoff relationships were detected as events in catchments with pervious lithology showed the highest correlation and the lowest scattering (Figure 3b). Rainfall events > 55 mm generated a R > 4 mm, excepting events occurred in late spring (Figure 3a) in catchments with predominance of pervious lithology (Figure 3b).

Event rainfall-runoff relationship was carried out considering the differences in the main land uses of studied catchments (Figure 3c). The median event R in forest (37.6 mm) was higher than agricultural (5.0 mm), shrub (2.0 mm) and agroforest (1.5 mm) catchments. However, the highest correlations between rainfall and R were obtained in catchments with a predominance of forest (R2 = 0.96), agricultural (R<sup>2</sup> = 0.61) and agroforest (R<sup>2</sup> = 0.58) land uses.

#### *3.2. Hydrological Response at Multiple Temporal Scales in Es Fangar Catchment*

#### 3.2.1. Annual Scale

The mean annual rainfall (840.8 mm yr−<sup>1</sup> <sup>±</sup> 213 mm y−1) calculated over the five hydrological years (2012–2017) was broadly representative (−9%) of the long term mean annual rainfall (927 mm y−<sup>1</sup> <sup>±</sup> 215 mm y<sup>−</sup>1: 1965–2016), but showed a high range of variation (coefficient of variation 25%, being 553 mm y−<sup>1</sup> to 1090 mm y<sup>−</sup>1), characteristic of Mediterranean conditions. Figure 4 shows the annual variability of rainfall, baseflow and quickflow contributions as well the annual Rc and the number of days with observed flow at the gauging station. Linear positive relationships were observed between annual rainfall, R, Rc, baseflow and quickflow (R2 <sup>≥</sup> 0.84, data not shown). The annual Rc ranged from 2.9% to 14.2% (mean value = 10.4%) and quickflow contribution from 9.9% to 45.0% (mean value = 33.0%). During the study period, flow was observed 42.8% of the time. From year to year, the cumulated number of days with flow ranged from 37.8% (138 days) to 65.2% (238 days) of the total days (Figure 4). An inverse relation between annual R and the number of days with flow was established. On the one hand, hydrological years with the largest R (2013–2014, 2014–2015 and 2016–2017) showed fewer days with flow and a lower baseflow contribution (<60%). In these years, 50% to 62% of the annual R was reached in 5 or fewer days (Table 2). In addition, days with more R were always during autumn and winter. On the other hand, the hydrological years with the lowest R (2012–2013 and 2015–2016) showed 208 and 238 days with flow. The baseflow contributions represented 70% and 90%, respectively, and 50% of the annual R was reached in 10 and 17 days, also respectively. The contribution of the 5 days with more R did not exceed 28% and 37%, respectively, and these days with more R were distributed among the four seasons.

**Figure 4.** Baseflow and quickflow contributions (mm) in the total flow for each hydrological year (October to September) at Es Fangar Creek. The annual Rc is depicted in %, as well as the total number of days (d) with recorded flow at the gauging station. The total annual rainfall is also illustrated in brackets.

#### 3.2.2. Seasonal Scale

Figure 5 shows the monthly mean values of rainfall, reference evapotranspiration and minimum, median and maximum R values. The mean monthly rainfall values were broadly comparable to those observed for the long-term period, except for October, when rainfall during the study period was much lower (56% lower than long-term rainfall).

The seasonal dynamics of rainfall and evapotranspiration controlled the R response. Characteristic wet (winter) and dry (summer) periods alternated throughout the year, separated by transition periods (last autumn and early spring). Autumn was the rainiest season, with low evapotranspiration and flow observed at the outlet for 52.7% of the time, on average. The autumn mean Rc was 9.1% and ranged from 4.1% to 11.4%. During November and December, a continuous baseflow was generated in response to a large rainfall amount after the pre-filling of the initial water storage in October. From October to December, the mean Rc increased from 1.3% to 13.1%. Finally, the minimum, median, average and maximum monthly R amounts in December were higher than in November, even with lower rainfall amounts (Figure 5). Although differences between the hydrological years were observed, autumn was the season with less inter-variability in terms of R response.

During winter, flow occurred during 90.6% of the time, although the mean rainfall amount (294 mm) was lower than in autumn (326 mm). This more frequent presence of flow was caused by the low evapotranspiration demand, thus keeping the hydrological pathways active, as already established in autumn. Consequently, in winter, a R amount similar to that of autumn can be generated

from a lower rainfall amount. The winter mean R coefficient was 16.9%, ranging from 1.5% to 27%. During winter, the mean monthly Rc decreased from 18.9% to 14.6% between January and March as a consequence of a decreasing rainfall amount and an increasing evapotranspiration demand.

In spring, flow was observed 41.3% of the time. The monthly rainfall decreased during the season and was much lower (127.4 mm) than in autumn and winter. Although monthly evapotranspiration losses increased and were higher than rainfall, the water reserves accumulated during the autumn and winter months sustained the flow contribution. The spring mean Rc was 5.1%, ranging from 0.4% to 8.4% during the study period. The average monthly R decreased during the season from April (7.7%) to June (0.3%).

Finally, in summer, flow was observed only 0.9% of the time. The stream was most often dry in this season due to the high negative balance between rainfall and evapotranspiration. R was only ephemeral in response to the episodic rainfall events more frequent in September, when 80% of rainfall occurred on average (Figure 5). The summer Rc was 1.4%.

**Figure 5.** Monthly time series of rainfall, R, reference evapotranspiration during the study period (2012–2017) at Es Fangar Creek. Box plots show minimum, median and maximum monthly R. Blue dots show mean monthly R. Long term (1964–2017) monthly rainfall distribution is also depicted.

#### 3.2.3. Non-Linearity Assessment at Event Scale

During the study period, the number of flood events per hydrological year was between 2 (2015–2016) and 14 (2014–2015). Their seasonal distribution was 13 events in autumn, 25 in winter, 7 in spring and 4 in summer. The events' characteristics are detailed in Table 3. An investigation of the variability and the non-linearity of the hydrological response at the event scale was carried out by means of a comparative assessment of seven rainfall-runoff events with similar Ptot, ranging from 41.8 to 49.8 mm, but different antecedent conditions and rainfall dynamics (Figure 6).

**Figure 6.** Selected events for non-linearity analysis at Es Fangar Creek. Events with a total precipitation between ca. 40 and 50 mm. The numbers located at the peak of each hydrograph indicate the ID of the events recorded in Es Fangar Creek during the study period 2012–2017, further exposed in Table 3.

The event on 21 January 2017 (Figure 6a) occurred under wet conditions following a runoff event on the day before. AP1d was 35.4 mm and Q0 was high (0.229 m<sup>3</sup> s−1). The rainfall characteristics (duration and intensity) were similar to those observed for other events of the same magnitude (Figure 6b,e,g), but the hydrological response was an order of magnitude higher in terms of R (28.9 mm), Rc (65%) and Qmax (3.942 m3 s<sup>−</sup>1) due to wet antecedent conditions (Table 4). This event had a similar duration to events that occurred under higher rainfall intensities (11 November 2012, 29 September 2014). However, wet antecedent conditions were much more favourable than rainfall intensity to reach a high Qmax and especially a higher R contribution.

The event on 4 April 2012 (Figure 6b) occurred in early spring, ergo at the beginning of the transition period from wet to drier conditions. The antecedent conditions were favourable to R generation, due to the presence of a baseflow maintained by water reserves accumulated in the autumn and winter months. Rc was 16%. The event on 11 November 2012 (Figure 6c) recorded the highest IPmax30 (i.e., 77.5 mm h<sup>−</sup>1) of the study period, occurring in the transition period from dry to wet conditions (middle autumn). Previous to this event, the stream channel was dry. Under these conditions of rainfall intensity and without previous baseflow, the Rc reached 14%. Although the event on 4 April 2012 recorded a similar Rc (16%) to the 11 November 2012 event, both events resulted in different hydrographs because of different rainfall characteristics (in terms of duration and intensity) and antecedent wetness conditions. As a result, the duration of 44.7 h for the event on 4 April 2012 contrasted with the duration of the 11th November 2012 event, which only lasted 18 h. In addition, the differences between these two events were also relevant in terms of Qmax—1.006 m3 s−<sup>1</sup> for the April event and 1.747 m<sup>3</sup> s−<sup>1</sup> for the event on 11 November 2012.

The event on 8 February 2017 (Figure 6d) occurred during the wet season. AP1d and a Q0 were, respectively, 0 mm and 0.012 m3 s<sup>−</sup>1. On the contrary, the event on 19 December 2016 (Figure 6e) occurred at the end of the transition period and had an AP1d and a Q0 of 2 mm and 0.012 m3 s−1, respectively. Both events had similar values of IPmax30, AP1d and Q0 (Table 4). Nevertheless, R and Rc were slightly larger for the event on 8 February 2017, likely due to the large rainfall amount accumulated during the two events (more than 400 mm of rainfall between 19 December 2016 and 8 February 2017).

The event on 7 May 2016 (Figure 6g) was similar to the event on 19 December 2016 (Figure 6e) in terms of Ptot, duration and intensity. Nevertheless, R response started earlier (i.e., shorter response time) for the event on 19 December 2016. This difference was caused because the December 2016 event occurred during the wetting-up period when catchment water reserves were increasing. As a result, whilst the R response was relatively small (4.4 mm), Rc (9%) and Qmax (0.394 m<sup>3</sup> s<sup>−</sup>1) were significantly larger than for the event on 7 May 2016 (Table 4), which occurred during the transition period towards dry conditions. With an AP1d of 1.6 mm and Q0 of 0.006 m<sup>3</sup> s<sup>−</sup>1, the lowest R response was generated (see total runoff, Rc and Qmax in Table 4). This R response was also the most delayed, as shown on Figure 6g, starting when 95% of the rain (39.8 mm) had fallen on the catchment.

The event on 29 September 2014 (Figure 6f) occurred during the dry season (i.e., late summer, early autumn) but with an AP1d of 11.8 mm. The rainfall event was the shortest of these seven selected events (3.2 h), with the second highest IPmax30 (56.4 mm h−1) of the study period. This very high intensity triggeredaRresponse during dry conditions. Similar rainfall but lower intensities during dry conditions always produced lower responses (Rc ≤ 3%). Although the R response was relatively small (4.2 mm) with Rc < 10%, Qmax was relatively high (2.356 m<sup>3</sup> s<sup>−</sup>1), promoted by the rainfall intensity.

#### **4. Discussion on Hydrological Responses in Small Mediterranean-Climate and Es Fangar Catchments**

#### *4.1. Annual and Seasonal Scales*

#### 4.1.1. Small Mediterranean-Climate Catchments: Lithology Influences

Catchments with annual rainfall ranging from 500 to 900 mm yr−<sup>1</sup> showed annual R ca. 20–350 mm, where the role played by lithology on the annual R response was important. Accordingly, the highest annual R values (>350 mm) were observed in catchments with badland areas [34] and a high subsoil

clay content promoting lateral flow and a perennial regime [35]. Swarowsky et al. [54] also reported an annual R 17% lower than in an adjacent catchment to that investigated by Lewis et al. [35] and related this difference to greater deep infiltration. Similarly to our results, Merheb et al. [46] obtained a significant rainfall-runoff correlation (R<sup>2</sup> = 0.69) using 160 catchments from the Mediterranean Region (0.35 to 21,700 km2). These authors highlighted that the scattering observed was due to karstic catchments or snowmelt contributions.

Also, scattering can be observed in the relationship for catchments with annual rainfall >1000 mm yr<sup>−</sup>1. For this group of catchments, the annual R ranged from 135 to 1200 mm, except for the Santa Magdalena catchment, which is entirely on a limestone lithology and yielded an annual R value of 70 mm [47]. Large annual R values usually corresponded to watertight catchments with R values > 1000 mm yr−<sup>1</sup> for two small Mediterranean-climate catchments on granite bedrock [36]. However, karstic environments can also promote very high annual R values, e.g., Ref. [55], as a result of incoming karstic spring sources. The conceptual model developed by Borg Galea et al. [2], which is based in exogenous and endogenous variables that link natural and human-derived variables and their influences in the system on a spatial–temporal scale, emphasized that climate is the main exogenous driver of the rainfall-runoff relationship, although their processes through precipitation, temperature and wind are mediated by catchment geology (endogenous variable).

Different seasonal rainfall-runoff relationships in relation to the influence of rainfall and evapotranspiration [56] are related to the alternation of rainfall and reference evapotranspiration throughout the year reproducing wet (winter), dry (summer) and transition periods (last autumn and early spring) [57]. Accordingly, winter is the season with the highest Rc, from 17% to 56%, due to rainfall being accumulated during autumn—and maintained in winter—and low evapotranspiration demand [10,15,56]. The same pattern was also observed during the five hydrological years assessed in Es Fangar Creek. In addition, all these studies also emphasized that in summer the hydrological response is limited or null, with Rc < 10%. Serrano-Muela [56] observed runoff coefficients < 5% in a range of rainfall amounts of 15–200 mm during summer in a small catchment characterized by a large forest cover. The driest environments in these studies obtained a summer Rc < 2%. Rc in autumn and spring were similar, ranging from 9% to 25% and 5% to 28% respectively. However, a variability throughout the seasons can be observed, especially in late autumn and early spring, as these are transition periods. Additionally, Lana-Renault [10] observed that the wetting-up period was steeper than the drying-down period, which was more progressive in a catchment located in the Pyrenees. Autumn was the season after the driest period, and it was when the evapotranspiration reached the maximum values. The beginning and the end of this season can be quite different in relation to Rc. The catchments have a null or limited response until a succession of rainfall events that fill the water reserves and generate favourable conditions for R generation [56]. The findings of Lana-Renault [10] in spring concluded that the accumulated rainfall during autumn and winter maintain high water reserves, allowing high runoff coefficients, even though spring was not the rainiest season.

#### 4.1.2. Es Fangar

In the Es Fangar catchment, a low mean annual R value (87 mm) was recorded during the 2012–2017 period, despite its mean annual rainfall value (black dot on Figure 2). This result is most likely related to the presence of massive karstic dolomites and breccias in the lithology of the catchment headwaters.

The reported annual R in Es Fangar is within the R range of small catchments located in Mallorca, with annual rainfall ranging from 500 to 900 mm yr<sup>−</sup>1, which yielded annual R between 36 and 130 mm. From these catchments, the highest R value was observed in an agricultural small lowland catchment prone to receiving subsurface contributions [15]. Spatial differences in the annual R between the headwater (102 mm) and catchment outlet (36 mm) due to transmission losses along the main channel were detected in the Sa Font de la Vila River, a mid-mountainous terraced catchment affected by forest fires [58]. Within this catchment, temporal differences were observed as Rc ranged from 1% to 22% from one year to year.

#### *4.2. Event Scale*

4.2.1. Seasonality, Lithology and Land Use Influences on Small Mediterranean-Climate Catchments Hydrological Response

The seasonal dynamics of rainfall and evapotranspiration during the year generate wet, dry and transition periods, which control the R response [10,53]. Therefore, seasonally, the highest correlation observed in spring could be related to the water reserves accumulated during the autumn and winter seasons, which is favourable for R generation [13]. According to the seasonal distribution of the 203 events, other studies in Mediterranean-climate catchments observed an event seasonality in the hydrological response [8,10,59] where autumn rainfall events produced a limited or null R response, especially at the beginning of the season. In addition, in winter and spring, the Rc values were always larger than 3%, being the maximum in spring (70%).

The role of transmission losses due to lithology has been investigated [32,60], depicting how low runoff values were recorded due to infiltration. In this study (see Section 3.2), events assessed in catchments with a predominance of pervious lithology showed a higher scattering than events in impervious catchments. In addition, rainfall events >55 mm generated R > 4 mm, except in events occuring in late spring in catchments with pervious lithology. In this case, the transmission losses due to pervious lithology have more influence on the hydrological response than the seasonal role of the water catchment reserves. Ries et al. [61] identified a low and large R response in rainfall events >50 mm based on the physiographic catchment and rainfall characteristics. Events occurring in catchments, which even had different land uses, with less bedrock permeability and less soil water storage but with higher values of rainfall intensity generate a greater R response than catchments with more bedrock permeability and soil water storage under low rainfall intensities. Contrarily, other types of lithology (i.e., badland areas) are characterized by a low infiltration capacity [62]. Latron and Gallart [63] observed that frequently badlands are the only active area to R response. According to these findings, Nadal-Romero et al. [64] observed an event stormflow coefficient of up to 20% in a Mediterranean-climate catchment with badland areas, even when the headwater catchment was forested (30% of the catchment area).

Events occurring in agricultural catchments obtained the second highest values in the median R and in the rainfall-runoff correlation, as arable land and its spatial distribution within a catchment have been demonstrated to be a driver for R generation [65]. However, events in the forest catchment had the highest correlation and median R, as the water yield and rainfall-runoff relationship decrease in forests catchments when compared to agricultural catchments, due to forest cover [66]. However, this correlation (R2: 0.96; n: 3) and median value is related to the presence of karstic springs in the forest catchment, which was the only forest catchment found in the literature. The median R of forest catchments could be similar to the values obtained in afforested ones, as both land uses trigger a reduction in the R generation. This reduction was quantified in a loss of water yield ca. 10–40% [38,66–69]. Nevertheless, larger flows and Qmax were observed in afforested catchments than forest ones. The shrub median R (2.0 mm) was higher than afforested catchments, as shrub coverage is characterized by lower vegetation cover than forest and afforested land uses. In addition, the second main land use in shrub catchments is agricultural land use (Figure 1A), which could explain the higher median R compared with the afforested catchments. The same results were observed in García-Ruiz et al. [66], where catchments under shrub land use had lower R contributions than catchments with forest land use due to a lower plant cover density in the shrub catchment.

#### 4.2.2. Rainfall-Runoff Relationship at the Es Fangar Catchment

The identification of driving factors explaining the variability of Rc can be explored through the relationships between Rc and R and Ptot, rainfall, maximum intensity and baseflow at the beginning of the event, as Latron et al. [8] carried out in a mountainous Mediterranean-climate catchment in North Eastern Spain.

Only the correlation between R and Rc was significant in the Es Fangar catchment (Figure 2a), explaining 73% of the variance. As in Latron et al. [8] this relationship had the highest correlation. However, contrarily to Latron et al. [8], the relationship between rainfall and Rc was highly non-lineal in Es Fangar, presenting a huge scattering, especially for Ptot between 20 and 60 mm (Figure 2b), which yielded Rc ranging from 1% to 65%. For larger rainfall events (i.e., 60 to 100 mm), the Rc values were high (>20%), but low R responses (with Rc ranging 1% to 3%) were also observed when large rainfall events were recorded after or during dry conditions. In line with the results from Latron et al. [8], Rc was not related with rainfall intensity (Figure 2c).

The relationship between Q0 and Rc was not linear in the Es Fangar catchment (Figure 2d). A threshold was identified leading to larger Rc (>10%) for events with Q0 above 0.04 m<sup>3</sup> s−1, which always occurred between November and April. The relationship between Q0 and Rc in Latron et al. [8] was stronger, because of the presence of a baseflow during most of the year. In general terms, the frequent absence of a baseflow in the Es Fangar catchment between May and October worsens most of the relationships presented in Figure 2 in comparison to a catchment with an almost permanent baseflow, as in Latron et al. [8] and Nadal-Romero et al. [40]. Then, baseflow–Rc correlation increased as the environment got wetter, from R2 = 0.19 in Es Fangar to R2 = 0.70 in a humid catchment [70].

In order to help in identifying the main factors involved in the hydrological response in the Es Fangar catchment, multiples relationships were carried out to investigate the influence of rainfall, Q0 and IPmax30 on R and Qmax (Figure 7). The hydrological response by using R and Qmax was highly non-linear, as observed in other Mediterranean-climate catchments [8,11] because the soil moisture antecedent conditions are crucial for R generation. Figure 7a shows how the R response could not be explained by Q0, as the highest R values presented a large scattering along the Y axis. This pattern confirmed that large R amounts (i.e., >10 mm) may occur from large rainfall events occurring in dry or wet conditions (i.e., with low or high Q0 values). Therefore, Ptot was a key factor for R response, as shown by the significant positive correlation (p < 0.01) between both variables (Table 2). The largest R values were always recorded during last autumn or winter, with contributions of >17 mm (Figure 3 and Table 5), and are characterized by Ptot, R, Rc, Qmax and AP1d values much higher than median values (Table 3). Contrarily, large rainfall events occurring in spring resulted in a small R contribution because of the dry antecedent conditions (i.e., low or null Q0 values). Other authors also assessed relationships between seasonality and R generation. In a mountainous catchment, Lana-Renault [10] observed that the highest rainfall-runoff correlations at the event scale were depicted in winter and spring due to higher water reserves than in autumn and summer. In Es Fangar, a group of events during autumn and winter characterized by Q0 > 0.080 m<sup>3</sup> s−<sup>1</sup> always generated Rc greater than 14% (Figure 2d). The Ptot of these events ranged from 5.2 to 61.4 mm but the AP1d was larger than 20 mm, except in one event (AP1d: 3.4 mm), which was characterised as the highest IPmax30 (i.e., 77.5 mm h<sup>−</sup>1) of the study period. These event characteristics suggest that Ptot was the main factor controlling an effective hydrological response as antecedent conditions (i.e., soil moisture degree of the catchment) have been demonstrated to play an important role in R generation [18]. Similarly, Lana-Renault [10] identified a cluster of events characterized with the highest Q0 that generated a larger R response, but the convective rainfall occurring in the Pyrenees during summer blurred the seasonal pattern described in the Es Fangar catchment.

Figure 7b shows that high Qmax values were observed in response to large rainfall amounts (above 50 mm) or for flood events with a Q0 value higher than 0.080 m<sup>3</sup> s−1. Both factors had a combined effect on the generation of high Qmax values as shown by the significant correlation observed between Ptot and Q0 and Qmax. The role of wet conditions over Qmax was also observed in other Mediterranean-climate catchments where AP1d correlated significantly with Qmax [15,17,71].

**Figure 7.** Multiples relationships between rainfall and R variables at the event scale at Es Fangar Creek: (**a**) rainfall-antecedent discharge-runoff; (**b**) rainfall–antecedent discharge-Qmax relationship; (**c**) IPmax30-rainfall-runoff and (**d**) IPmax30-rainfall-Qmax. The points with an ID are depicting the selected events for the non-linearity assessment, as it can be further consulted in Figure 6 and Table 3.

**Table 2.** Pearson correlation matrix between selected variables.


R and Qmax were not significantly correlated with IPmax30 (Table 2). Figure 7c,d show how high values of R and Qmax resulted for events with low IPmax30 when rainfall was > 50 mm or high Q0. This absence of correlation between rainfall and the catchment response (R and Qmax) indicates that R generation is most likely the result of saturation excess processes [72], where Ptot and antecedent conditions played a key role. For events between June and September, saturation excess processes were, however, probably less dominant and combined with infiltration excess processes in response to

highest IPmax30 values ranging from 19 to 56.4 mm h<sup>−</sup>1. As a result, the events observed between June and September were shorter (5.7 to 14 h duration compared to 8 to 65.5 h for wet conditions events). Therefore, R mechanisms co-exist within a Mediterranean catchment [13] such as Es Fangar, with a marked seasonality, with saturation processes being dominant during the wet period and Hortonian R during the dry period. Consequently, the alternation of the R mechanisms generated lower and larger R contributions in events occurring during the dry and wet seasons, respectively. These processes caused the non-linearity of the hydrological response because the same amount of Ptot generated a non-lineal R response due to different soil moisture conditions [15,64].

#### 4.2.3. Non-Linearity Assessment at Es Fangar Catchment

The events analysed in Figure 6 presented a huge variability in their hydrological response, despite recording a similar Ptot, ranging from 41.8 to 49.8 mm. Rc ranged from 1% to 65%, depending on the catchment moisture conditions, rainfall intensities and seasonality characteristics. The highest hydrological response occurred under marked wet soil moisture conditions in the winter period, even with low rainfall intensities (Figure 6a). A similar complexity in R responses has already been observed in other Mediterranean-climate catchments where seasonality and antecedent soil moisture conditions played a key role in R generation [11,57]. Gallart et al. [73] also showed that different R generation behaviours could be observed during the year due to varying catchment wetness conditions and changing rainfall events characteristics. The results obtained by these authors are clearly illustrated by comparing Events 2 and 20 with Events 18 and 37 in Figure 6. Indeed, R events that occurred during wet conditions with low rainfall intensities (18 and 37) showed similar Rc than events that occurred during dry conditions with high rainfall intensities (2 and 20). However, during dry conditions, runoff events resulting from low intensity rainfall yielded the lowest R responses (e.g., Event 36). Our results are also in agreement with the observations made by Schnabel et al. [74], who observed higher R contributions for low rainfall intensities events occurring in periods with high soil water content than for events with high rainfall intensities. Besides, those authors pointed out that this condition was common during years with above average rainfall (i.e., during wet years) but was rarely observed during dry years.

#### **5. Conclusions**

The evaluation of multiple temporal scales in contrasting small Mediterranean-climate catchments has improved the understanding of the role played by lithology and land use in the hydrological response. The assessment of rainfall-runoff relationships at the annual scale in small Mediterranean-climate catchments showed a significant linearity in the hydrological response due to the importance of the annual rainfall amount. Nevertheless, lithology effects on R generation explained an increase in the scattering in the rainfall-runoff relationship because pervious and impervious materials triggered larger and lower R contribution, respectively. Despite the significant correlation between rainfall and R, the Es Fangar Creek dataset illustrated a huge intra-annual variability of the rainfall-runoff relationship during the five hydrological years analysed, as seasonal rainfall and evapotranspiration dynamics controlled the R response. These dynamics were clearly observed in the average seasonal Rc, decreasing from winter to summer. As a result, these seasonal differences should be considered as a starting point of the non-linearity generation in the rainfall-runoff relationships at the event scale.

At the event scale, lineal and non-lineal relationships were observed in the rainfall-runoff relationships in small Mediterranean-climate catchments, suggesting that different factors conditioned the R response. The total rainfall was the most significant driving factor, although the interaction between seasonality and the spatial diversity of lithology and land uses at the catchment scale also played an important role in R generation. Thus, the highest correlations at the seasonal scale were observed in those events which occurred in winter and spring when the highest water reserves favoured the R response. Lithology caused a higher dispersion in the rainfall-runoff relationships at the event scale in the set of small Mediterranean catchments because pervious materials required higher antecedent wetness conditions. Different land uses promoted a decrease in R generation, comparing agricultural with scrubland uses, because agriculture promoted the highest correlation in the rainfall-runoff relationships due to lower vegetation cover. Therefore, human-induced alterations related to land uses in the R generation must be assessed in those countries were the abandonment of upland agricultural catchments have allowed natural vegetation to expand (i.e., southern Europe) and also in countries were agriculture and irrigated areas are increasingly replacing forest and shrub lands (i.e., North African Maghreb). However, events under agricultural land use or which occurred in the winter season independently of the land use and lithology catchments were the situations which were able to generate a higher linearity in the rainfall-runoff relationships.

This temporal downscaling from the annual to event scale elucidated how different R mechanisms can co-exist in small Mediterranean-climate catchments, considering the main temporal and spatial factors that govern the river catchment connectivity. Despite this, controls on R generation in Mediterranean-climate catchments require further attention to assess the role of lithology, land use and seasonality and their combined effects on the hydrological response for going beyond in the comprehension of highly sensitive areas to global change, such as the Mediterranean region. Considering that the performance gaps of hydrological processes in hydrological models due to the results of rainfall-runoff at daily scales are smoothed in small Mediterranean-climate catchments, the findings of this study improve the comprehension of hydrological processes as practical field knowledge needed to complete hydrological modelling at the event scale.

**Author Contributions:** Conceptualization, J.F., J.L. and J.E.; methodology J.E. and M.T.-B.; data curation, J.G.-C., A.C. and M.T.-B.; writing-original draft preparation J.F.; writing-review and editing, J.F., J.L., J.C. and J.E.; visualization J.G.-C., J.C. and A.C. and supervision J.L. and J.E. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was supported by the research project CGL2017-88200-R "Functional hydrological and sediment connectivity at Mediterranean catchments: global change scenarios –MEDhyCON2" funded by the Spanish Ministry of Science, Innovation and Universities, the Spanish Agency of Research (AEI) and the European Regional Development Funds (ERDF). The contribution of Jérôme Latron was supported by the research project PCIN-2017-061/AEI also funded by the Spanish Government. Josep Fortesa has a contract funded by the Ministry of Innovation, Research and Tourism of the Autonomous Government of the Balearic Islands (FPI/2048/2017). Julián García-Comendador is in receipt of a pre-doctoral contract (FPU15/05239) funded by the Spanish Ministry of Education, Culture and Sport. Miquel Tomàs-Burguera acknowledges the support from the project CGL2017-83866-C3-3-R financed by the European Regional Development Funds (ERDF) and the Spanish Ministry of Science, Innovation and Universities. Jaume Company is in receipt of Young Qualified Program fund by Employment Service of the Balearic Islands and European Social Fund (SJ-QSP 48/19). Aleix Calsamiglia acknowledges the support from the Spanish Ministry of Science, Innovation and Universities through a pre-doctoral contract (BES-2013-062887).

**Acknowledgments:** Part of the meteorological data were provided by the Spanish Meteorological Agency (AEMET).

**Conflicts of Interest:** The authors declare no conflict of interest.


*Water* **2020**, *12*, 299

**Appendix** 

**A**


**Table A1.** *Cont*.


**Table 2.** Relative rainfall and R contribution of the highest 5 days and number of days to reach the 50% of R for each hydrological year in the Es Fangar Creek.


**Table 3.** Flood event characteristics during the study period 2012–2017 in the Es Fangar Creek.


**Table 4.** Main characteristics of the selected events for non-linearity analysis in Es Fangar Creek.

**Table 5.** Main characteristics of the highest runoff event contribution in Es Fangar Creek during the study period 2012–2017.


**Figure 1.** Land uses (%) for the 10 of 12 small Mediterranean-climate catchments selected for the rainfall-runoff relationship analysis. Guadalperalón and Parapuños catchments have not detailed land uses information.

**Figure 2.** Relationship between (**a**) R and Rc coefficient, (**b**) rainfall and Rc, (**c**) IPmax30 intensity and Rc and (**d**) base-flow specific discharge and Rc at Es Fangar Creek. Dotted lines show significant (p < 0.01) fits with a power function.

**Figure 3.** Highest R event contribution during the study period (2012–17) at Es Fangar Creek. The numbers located at the peak of each hydrograph indicate the ID of the events recorded in Es Fangar Creek during the study period 2012–2017, further exposed in Table 3.

#### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Water* Editorial Office E-mail: water@mdpi.com www.mdpi.com/journal/water

MDPI St. Alban-Anlage 66 4052 Basel Switzerland

Tel: +41 61 683 77 34 Fax: +41 61 302 89 18