Modeling and Upscaling Plot-Scale Soil Erosion under Mediterranean Climate Variability

Abstract

Soil erosion is an issue in the Mediterranean slopes. Erosion plots are useful to quantify erosion rates, but data are difficult to scale up to a slope level. Moreover, short observational frameworks are generally established, making it difficult to represent multi-year fluctuations. This paper deals with the potential of parsimonious modelling to upscale plot erosion (~23 m²) at Monte Pino Met European Research Observatory (South Italy) from 2001 to 2006. Under the assumption that the slope is fractal and contains plots, monthly gross soil erosion was modeled by lumping together the erosivity factor (runoff component), Normalized Difference Vegetation Index (vegetation cover factor), and the spatial scale dependence (slope length factor). This model was applied to reconstruct monthly gross soil erosion rates for the period of 1986–2006, for which hydrological inputs were available with sufficient detail. Pronounced interannual variations, with two distinct patterns, were observed: increasing rates of erosion were visible in 1995–2006 (peaking in November 1997, 50 Mg·ha⁻¹·month⁻¹), while in previous years only a few peaks slightly exceeded the average of the whole period (1 Mg·ha⁻¹·month⁻¹). Hydrological conditions indicate that important erosional processes have been triggered during low-frequency, short rainfall events occurring in spring–summer (e.g., May 2001, June 2003), or during longer, less intense events occurring in autumn–winter (e.g., November 1997) seasons. It is likely that increased precipitation amounts associated with more frequent convective storms created conditions for higher energy events triggering erosion. For the recent warm period, investigations at a higher than monthly resolution are required to better assess the seasonal changes of erosion rates and their relationship with soil conservation.

1. Introduction

Soil erosion is a form of land degradation caused by complex and not always well known interactions among rainfall impact and variability, overland flow, land-use and land-cover, soil properties, topography, and conservation management [1]. In agricultural areas, tillage impact is noticeable on soil’s physical, chemical, and biological properties, which includes the effect on the soil environment in the form of runoff and soil erosion [2,3]. Extreme precipitation events can lead to considerable soil loss, with implications for downstream flood risk, and shallow landslides [4]. Soil erosion is a great concern all over the world, and a primary problem in many Mediterranean-climate areas of Southern Europe, which are characterized by typical within-year seasonality in precipitation distribution, with heavy rainfalls in fall and spring, and a relatively large year-to-year climate variability [5,6]. This feature affects, in turn, the local hydrological conditions with related instabilities in the vegetation ecosystems when severe perturbations occur [7]. In Mediterranean environments, soil erosion is mainly driven by raindrop impact and overland flow, and induces a loss of cropland available for food production, and is a major cause of various forms of land degradation [8,9,10,11]. Such meteorological and hydrological events exacerbate soil erosion phenomena, especially in those locations where the soil parent material is formed by a weakly consolidated substrate, and observed changes in local climate conditions reveal an increase in rainfall intensity [12,13,14]. A rapid evolution towards the conditions of severe soil erosion degradation is ongoing in the Mediterranean belt [15,16,17], including Southern Italy [18,19,20].

Soil erosion assessment and mitigation requires large datasets to understand the fluctuations in time and space, and is particularly challenging under climate change-driven shifts in precipitation [21]. For the Euro-Mediterranean region, only long-term averages of sediment yield have been predicted by the PESERA (Pan-European Soil Erosion Risk Assessment) project (see Figure 1; [22]). Available studies about the recurrence of dry periods and extreme rainfall events, which affect the variations in monthly soil erosion rates [23], should be extended over longer time intervals, and carried out at different temporal scales (e.g., inter-annual and inter-monthly scales) in order to be able to highlight not only the differences between the spring–summer and autumn precipitation, but also the seasonal variability in rainfall frequency and intensity. Month-time intervals seem an appropriate timescale to identify both short-and-long-term changes in soil erosion and its effects on ecosystem evolution. However, the complex interaction that exists between factors governing soil erosion processes (e.g., soil physical and hydraulic properties, local terrain attributes, and land-use) makes soil erosion predictions very challenging [24]. Whereas the available literature and databases are sources of good climatic data at regional and sub-regional scales for the last century, they do not provide enough information on space-time scale, location, and timing of geomorphologic changes (except in the general sense of rivers or coastline modifications [25]). Soil losses were quantified in various Mediterranean landscapes by using erosion plots [26,27,28,29,30,31]. The results of these studies have provided important reference outcomes that have had a fruitful repercussion in a range of related disciplines, such as ecology, geomorphology, and hydrology [32,33]. However, monitoring of soil erosion phenomena is often restricted to a few years of observations only, implying that the measurements are not necessarily representative of middle-term fluctuations occurring in ecosystems, especially under a Mediterranean climate [34,35].

We investigated the possibility of reconstructing middle-term monthly gross soil erosion on the basis of discrete plot measurements and subsequent upscaling of plot-scale modeling results to a slope-scale. This is a key issue for soil loss predictions also considering that soil erosion, as well as other environmental phenomena, is often monitored only across restricted areas. For this purpose, an empirical model is developed and calibrated at a plot-scale, and then evaluated to assess its effectiveness of transfer information to the slope scale though a GIS-integrated procedure. The model is used to reconstruct monthly gross soil erosion rates for the period of 1986–2006, and investigates possible signals of monthly erosion change at the Monte Pino plot near the city of Benevento in Southern Italy (Figure 2). Climate extremes are also considered. Under the conditions of climate change, a divergence of precipitation patterns is predicted. Some regions of Mediterranean Central Area (MCA) will experience more annual precipitation and a more regular seasonal distribution of rainfall events, while other areas will experience stationarity in annual amounts of precipitation and a higher rainfall intensity.

2. The Monte Pino Experimental Site (MPES)

The Monte Pino experimental site (MPES) is located in the Calore River Basin, in the Samnium territory (14°42′30′′ long. East and 41°06′30′′ lat. North, 184 m above the sea level; Benevento, southern Italy; Figure 2a,b). It is placed between the Benevento valley (northeast) and Montemauro Mountain (southwest; dotted square in Figure 2c) along the eastern slope of the Monte Pino (hence the name of the experimental site). The MPES comprises of an experimental plot of 23 m² and a meteorological station. The meteorological station is a part of the Met European Research Observatory (MetEROBS) of the HyMex network [37], and is located, together with the erosion plot, in the extreme east mountainside (little white square in Figure 2c) that is also part of the regional park of the Taburno and Camposauro Mountains. The HyMex project develops and facilitates hydrological monitoring in European areas. This site offers a unique opportunity for modeling erosion responses to climatic variability in an area that has not been subject to drastic land-use change over the last century.

The images of Figure 3 illustrate the eastern slope of MPES as a heterogeneously structured entity, having a common climate (Figure 3a) with complex pedo-topographic and geological features, and a variable land-use (Figure 3b) controlling their hydrologic dynamics. The area has a typical Mediterranean climate, which allows growing high-value crops such as grapes and olives to produce high-quality wine and olive oil. In particular, mean annual rainfall is 950 mm (±150 mm Std.Dev.), and mean annual temperature is 14 °C (±0.5 °C Std.Dev.), both computed using the available datasets from 1986 to 2006. The MP site has a surface area of about 10 km² and a mean elevation of 275 m above sea level (Figure 4a). The mean slope steepness is 16%, with hillsides mainly oriented towards the northwest, sloping to around 30% (Figure 4b). A plot of about 23 m² was established along the hillslope, in an area characterized by a slope gradient of 32% and cultivated with olives trees (Figure 4c).

3. Data and Model

3.1. Data

We used climatic data retrieved from 1986 to 2006 from the meteorological station, installed within the Monte Pino experimental site, and soil erosion rates collected from 2001 to 2006 at the 23 m² experimental plot. Rainfall intensity was measured using a tipping bucket rain gauge every 24 h, whereas the soil erosion rate was measured after each rainfall event. The tipping bucket was installed on a pole that was 1 m from the ground. The plot was tilled once each year (usually in spring) to reproduce the MPES land-use. This also allows for reducing the risk for the plot to exhaust the available material for soil detachment over long-term studies. Runoff was drained into a sedimentation drainage installed in the ground in the lower part of the plot. Eroded soil (kg·m⁻²) was collected once per month. Measured soil erosion rates are reported in

Table 1. Summary of rainfall and erosion rates for each rainstorm event between May 2001 and November 2006 recorded at Monte Pino—MetEROBS.

Year	Month	Day	P (mm·day−1)	EI30 (MJ·mm·ha−1·h−1)	Erosion (kg·m−2)	Events
2001	May	23	67.0	2300	20.000	Downpours with flash-floods and surface landslides
		24	26.6	200	1.740	Thunderstorm
	August	22	17.0	177	0.668	Thunderstorm with isolated downpours
	September	15	14.8	33	0.056	Thunderstorms and showers
		17	10.0	50	0.160	Intense rainfall
	November	13	50.0	261	0.800	Intense rainfall with wind gust
		28	18.0	40	0.223	Rainfall
2002	April	4	41.0	47	0.044	Intense rainfall
	May	20	13.0	88	0.018	Thunderstorms and showers
		27	8.6	47	0.031	Short thunderstorm
	June	2	17.2	65	0.013	Thunderstorms and showers
	September	22	66.0	494	0.578	Continuous thunderstorms and showers with wind
		23	26.6	113	0.036	Thunderstorms and showers (afternoon)
	October	10	46.4	171	0.033	Showers at morning and thunderstorm at evening
2003	April	6	26.6	20	0.016	Continuous and moderate rain, somewhat with showers
	June	3	60.0	1500	2.990	Afternoon deluge with quiet rain, after stormy
	July	7	17.0	80	0.044	Moderate rain with isolated thunderstorms in the afternoon
		31	27.2	172	0.333	Heavy thunderstorm at evening
	October	5	29.6	188	0.261	Heavy thunderstorm followed by moderate rain at nighttime
		29	27.8	136	0.036	Thunderstorms in the evening
2005	August	21	34.0	86	0.267	Afternoon thunderstorm
	September	1	14.0	80	0.220	Afternoon thunderstorm
		18–20	32.7	60	0.600	Showers
2006	July	25	18.2	186	0.187	Evening thunderstorm
	August	6	10.8	38	0.049	Afternoon thunderstorm
		8	20.0	148	0.190	Afternoon thunderstorm
		13	38.2	324	1.320	Afternoon thunderstorm
	September	15	34.8	218	0.222	Showers
		26	39.6	277	0.266	Thunderstorm and showers
	November	22	32.5	45	0.022	Showers and isolated thunderstorms

P: Precipitation; EI₃₀: Rainfall erosivity.

. Both sediment and rainfall data were aggregated on a monthly basis. For instance, soil loss measured in the plot, during the years 2001 and 2002 was of 28 kg·m⁻² and 0.55 kg·m⁻², respectively. From 2001–2006 the median erosion values were equal to 2.52 kg·m⁻². This example indicates the importance of sparse disastrously rainy days or months in the determination of medium-term erosion, and the need to explore the entire distribution of the erosion rates by months over a long period. Model upscaling was completed into a GIS environment using data derived from a 20 m² square-cell digital elevation model (DEM) and a color orthophoto. The DEM was used to derive topographic attributes (Figure 4b), whereas the orthophoto was used to obtain values for the Normalized Difference Vegetation Index (−1 ≤ NDVI ≤ 1), a metric of vegetation greenness (usually derived from satellite-driven information), with higher values indicating highly vegetated areas.

3.2. Model

The Universal Soil Loss Equation (USLE) and its modified versions, such as RUSLE or MUSLE, are definitely the most known and used regression-based soil erosion models [38,39,40], and also provide the basis of process-based hillslope models, such as the USDA Water Erosion Prediction Project (WEPP; [41] and Eurosem [42]). Process-based models are data-demanding, and their applications at relatively large spatial scales can become a drawback developed a model [43] that is aimed at predicting monthly soil erosion caused by rainfall and overland flow across a plot-cell of 25 m × 25 m and is cast as follows:

$$E=k\cdot Q^2\cdot S^{1.67}\cdot e^{-0.07\cdot v}$$

(1)

where E is monthly soil erosion (mm·month⁻¹), k is the erodibility index, Q is overland flow (mm·month⁻¹), S is slope steepness (m·m⁻¹), and v is vegetation cover (%).

In our study, a scale transfer tool is required to expand Thornes’s model [43] from the plot scale (23 m²) to the slope scale (an area of about 10 km²). Since the spatial averaging is basically equivalent to the spatial integration [44,45], a viable possible upscaling technique is by averaging the plot-scale model over the entire mountainside area [46]. This approach is employed in our study and based on the assumption that the plot is a fractal-element of the slope. Accordingly, Thornes’s model [43] is revised to generate an adapted RUSLE-Thornes model, which replaces runoff, Q, with RUSLE erosivity-factor, EI, and vegetation cover, v, with the Normalized Difference Vegetation Index, NDVI_m. Therefore, the time-scale dependence of erodibility was expanded to include spatial scale-dependence of slope length factor (L).

These attributes are lumped together into the following semi-empirical, nonlinear relation:

$$E_{\left(\mathrm{ART}\right)m}=\left[k_m\cdot {\left(E{I}_m\right)}^{\mathsf{\eta }}\cdot L\cdot \left(\gamma \cdot d_{om{s}^{1-D}}\right)\right]\cdot e^{\left(-\omega \cdot 100\cdot NDV{I}_m\right)}$$

(2)

where E_(ART)m is the gross soil erosion rate in kg·m⁻²·month⁻¹. In Equation (2), the roughness erosive-resistance (the exponential term in the right-hand side), described as a function of vegetation cover, is balanced by the interaction between the potential transport capacity (the term in square brackets), described as a function of soil, rainfall, and topography. For the sake of simplification, Equation (2) does not explicitly account for rainfall-runoff processes (with the so called rainfall-runoff erosivity factor, R, calculated as a summation of event-based energy-intensity values, EI₃₀), under consideration that in the USLE some relationship does exist between runoff and EI₃₀, and between the peak rate of runoff rate and the maximum rainfall intensity in 30 min (I₃₀) [47]. The variables between square brackets are explained as follows:

$$k_m=K\cdot \sqrt{W\cdot \left(d_{wet}+C\right)+Z}$$

(3)

where k_m is the erodibility index, d_wet is the number of rainy days in a month with rain depth ≥1 mm·d⁻¹, K, W, C, and Z are empirical parameters;

$$E{I}_m=0.1174\cdot \left(\sqrt{p_r}\cdot d_r^{0.53}\cdot h_r^{1.18}\right)$$

(4)

where EI_m is the rainfall erosivity-factor (MJ·mm·ha⁻¹·h⁻¹·month⁻¹) that provides the forces applied to the soil and causing water erosion (indirect RUSLE procedure; details in [48]). In Equation (4), p_r is the monthly precipitation amount (mm), d_r is the daily maximum rain in each month (mm), and h_r is the hourly maximum rain in each month (mm);

$$\left(\gamma \cdot d_{oms}^{1-D}\right)$$

(5)

is the variogram equation to calculate the fractal dimension of topography [49], which can be converted to the equality;

$$\left|\frac{Z_p-Z_q}{d_{oms}}\right|=\gamma \cdot d_{om{s}^{1-D}}$$

(6)

where Z_p and Z_q are the elevations at the locations p and q, respectively. d_oms is the distance between p and q (that can be viewed as the original modelling scale), γ is a coefficient, and D is the fractal dimension. The following equation:

$$\left|\frac{Z_p-Z_q}{d_{oms}}\right|$$

(7)

represents the surface slope, so it can be assumed that the slope value, S (in Equation (1)), is associated to its corresponding scale–grid size d_oms [49]. As the USLE was originally developed for agricultural fields, its application to landscape scale erosion modeling is often inappropriate. Developed models based on the unit stream power theory [50,51] include the influence of terrain shape and are more suitable for complex topographic conditions, permitting to accommodate issues of upscaling with soil erosion models. USLE and RUSLE, which consider soil erosion along the flow line (such as in a plot field), have disadvantages because they neglect the influence of flow convergence/divergence, and their application is limited only to areas experiencing net erosion. The general form of the sediment transport equation can be used to describe the effects of terrain and scaling on soil erosion [52]. In particular, the following term in Equation (2):

$$L\cdot \left(\gamma \cdot d_{oms}^{1-D}\right)$$

(8)

can be viewed as the RUSLE-slope length factor LS (m). According to Mitasova et al. [51], the unit stream power varies with the runoff contributing area (A, in m²·m⁻¹), and the slope angle β (degrees), as:

$$L\cdot \left(\gamma \cdot d_{om{s}^{1-D}}\right)=\left(m+1\right)\cdot {\left(\frac{A}{a_o}\right)}^m\cdot {\left(\frac{\sin \beta }{b_o}\right)}^n$$

(9)

where $\left(\gamma \cdot d_{oms}^{1-D}\right)$ is the variogram equation to calculate the fractal dimension of topography in the upscaling equation, a_o = 22.13 is the RUSLE unit plot length (m), b_o = 0.0896, and m is the slope length exponent equal to (β/(1 + β)) being related to the ratio β of rill erosion (caused by overland flow) to interrill erosion (principally caused by raindrop impact; [51]). Parameters m and n are set equal to 0.5 and 1.3, respectively [51]

The upslope contributing area is calculated from the Digital Elevation Model having a cell size of 20 m, and using the SINMAP tool of ArcView GIS [53,54]. This tool employs the D-∞ technique proposed by Tarboton [55], that computes the contributing area more accurately on divergent hillslopes and uses a multiple-flow direction algorithm allowing for (i) continuous flow angles, and (ii) flow partitioning between one or two neighbor pixels. For the plot-scale based model:

• the contributing area is equal to the plot area; and,
• the term e^{(−ω∙100∙NDVIm)} is the Thornes vegetation exponential function, where NDVI_m is the Normalized Difference Vegetation Index derived from the red near-infrared reflectance ratio (NDVI = (NIR − RED)/(NIR + RED)) (where NIR and RED are the amounts of near-infrared and red light, respectively, reflected by the vegetation and captured by the sensor of the satellite).

NDVI values were estimated indirectly according to [56]. NDVI takes on values of 0 for bare soil, and 0.8 for maximum vegetation greenness. Parameter ω is a function of the ratio β defined above [43]. Exponent η of Equation (2) depends on the type of sediments transported by the erosion process. In similar problems, Carson and Kirkby [57] showed that various suspended sediment and the surface erosion equations correspond to η = 1.0, 1.5, 1.6, and 2.0.

3.3. Model Calibration

The model for predicting monthly gross soil erosion rates was run at MetEROBS over the years from 1986 to 2006. The data gathered over the overlapping period between the modeled time series and the plot measurement (i.e., from 2001 to 2006) have been used to calibrate the model parameters. Model calibration was performed using the measurements collected in the 23 m² experimental plot of the site. Instead, model validation was carried out by comparing the modeled data with the datasets collected by the Italian Hydrographic and Geographic Service within the Calore River Basin (having a drainage area of about 3000 km²), during the years of 1959 and 1961 [58]. We used only these two years for model validation because of the absence of longer overlapping time-series.

The values of the empirical parameters were determined at the calibration step using mean-monthly measurements through a least-square minimization technique. An iterative calibration process was employed to identify the terms of Equation (2). First, the set of parameters was determined for k_m Equation (3), fitting the erosion values and keeping constant the parameters of Equation (2). Next, the parameters of Equation (2) were calibrated against the variables measured at the plot scale. The process was reiterated up to reach a converging solution. The parameterization obtained (see

Table 2. Calibrated parameters for the study area.

Parameter	Values
η	1.66
K	0.001
W	0.16
C	5
Z	1
ω	0.05

) roughly matched the soil erosion data of both calibration and validation datasets (Figure 5a,b). A set of performance statistics—Percent Relative Root Mean Square Error (RRMSE%, from 0, best, to positive infinite [59], Modelling Efficiency (EF), from negative infinite to 1, best [60], Coefficient of Residual Mass (CRM), from negative infinite, over-estimation, to positive infinite, under-estimation [61])—computed on the basis of the calibration dataset indicated a good agreement between modeled and measured data (RRMSE = 12%; EF = 0.81; CRM = 0.15).

4. Results and Discussion

Equation (2) is used to reconstruct the historical time variations of the monthly gross soil erosion rates along the eastern slope of the Monte Pino site. The bar graph in Figure 6 displays the amount of particles mobilized within the Monte Pino site, whereas the blue lines show the erosivity density as computed by the ratio of simulated erosivity, versus the observed rainfall depth from January 1986 to December 2006. In the graph showing the gross soil erosion time series, one can easily detect two main periods. The first period spans from 1986 to 1994, when the monthly gross erosion values take on values near the median and attain the maximum values close to the average. The second period spans from 1995 to 2006, and is characterized by extreme values of monthly gross erosion abruptly increased with the highest peak located in November 1997. Although such increments match the general trend of variation of the erosivity density for the same period, the November 1997 peak does not correspond to a peak of the erosivity density. This can be interpreted by the fact that during this precipitation event the larger amount of soil particles are mobilized by overland flow.

Gross soil erosion time series is depicted in Figure 6, and characterized by a non-Normal probability distribution with a relatively great standard deviation value (not shown). The difference between the mean and median values of this series indicates also an erratic distribution of erosion-events that, in particular months, might also exceed the tolerable soil loss long-term annual-threshold (about 7 Mg·ha⁻¹; [58]). This was the case of the events that occurred in January, August, and November 1997, May 2000, and 2001, June and September 2003.

It is important to note that our model runs at a monthly timescale and is not able to predict soil erosion rates caused by daily storm events. For instance, more than 70% of sediment produced in November 1997 derived from an erosion process that occurred on 13 November, when the rainfall depth recorded at MetEROBS was 111 mm, of which about 86 mm fell in only three hours. Due to this event, remarkable soil erosion, floods, and also landslides occurred in the entire Benevento Province. During this event, in the Serretelle River Basin (including the MPES) the upper parts of high order channels were extensively incised (laterally and vertically) in response to heavy rainfall (Figure 7a). For this area, the analysis of historical precipitation data, referring to three-hour rainfall (Gumbel method, [62]), shows that a rain event of 86 mm has a return period of 20 years. On 13 November 1997 the rainfall erosivity was 933 MJ·mm·ha⁻¹·h⁻¹, a value comparable with the annual average value recorded for this area. This gives an idea of the high erosive-and overland-flow potential of this storm.

Most of the intense rainstorms with a low spatial extent, including the event of May 2001, occurred between the end of spring and summer. In May 2001, a heavy downpour hit the agricultural lands of the MPES, inducing extreme topographical changes connected with triggered erosion-processes (Figure 7b). In particular, precipitation characterizing this meteorological event amounted to 65 mm in 35 min, with a maximum intensity in 30 min higher than 100 mm·h⁻¹. Field observations revealed important changes in the landscape caused by this event, like deepening of the rill network and development of new rills even in low and flat places. Since extreme rainfall events (e.g., greater than 100 mm·d⁻¹) are not frequent at MetEROBS, the soil material forming the hillslope surface is easily transported along the slope also as consequence of moderate events (rainfall intensity roughly less than 20 mm·h⁻¹), producing changes and exporting the material out of agricultural fields.

Figure 7 illustrates two situations in which remarkable land degradation phenomena at MPES have been caused by extreme rainfall events. On November 1997 and May 2001, which can be considered as representative months with respect to modeling erosion processes under Mediterranean climate conditions, the estimated sediment was close to 50 Mg·ha⁻¹ and 30 Mg·ha⁻¹, respectively. Such values are much more above the long-term monthly average of gross soil erosion of 4 Mg·ha⁻¹ and 2 Mg·ha⁻¹. Based on both field observations and modeling outcomes (i.e., see Figure 7, left) we put forward the hypothesis that during the event occurring on November 1997, soil particles moving within MPES were successively transferred to the floodplain via drainage channels, and finally into Serretelle river. Field observations, on the other hand, indicated that during the extreme event that occurred on 14 May 2001, the soil material mobilized along the Monte Pino slope was prevalently trapped within concave subsections of the area. Our interpretation is that the hydrologic condition of the upper pedogenized layer (i.e., water content), forming the ground surface, influence the persistence and duration of runoff (together with morphology and land-use/cover of the slope). In other words, during the November 1997 event runoff duration was longer than the time of concentration (t_c) for this part of the basin, and thus mobilized sediment that reached the river. Conversely, during the May 2001 event, runoff was concentrated within the MPES area, and its duration was shorter than the time of concentration. Thus, most of the mobilized sediment did not reach the river outlet.

On these bases, erosional processes in the MPES can be triggered by meteorological events of magnitude around the average (for the site), but also during low-frequency, high-magnitude short events. This is the case of the extreme meteorological events occurring in May 2001 and June 2003. Field observations indicated that during high-magnitude phenomena, most of the soil lost along the hillslope was deposited before reaching the main stream outlet. Therefore, soil erosion within the MPES was controlled by spring-summer high-intensity rainfall events, and autumn and winter prolonged rainfall inducing persistent runoff.

The seasonal variability in erosivity and gross erosion is very important for evaluating the seasonal hazard of erosive rainstorm events. Figure 8 highlights the seasonal variability that affects the ratio of monthly to annual erosivity and gross cumulative erosion. The relative proportion between the monthly erosivity and erosion ratios remains similar, except for April, August, September, November and December. Especially the months of September and November present a higher disproportionality connected to the ability of heavy rainfall events to trigger erosive processes, and to soil condition and soil use in these months.

Seasonal Timing of the Intensifying Storminess

A more enlarged view of the Mediterranean Central Area (MCA) over the recent decades shows that autumn seasons are prone to intensified precipitation rates, with positive anomalies over many zones of central and southern Italy. Increased temperatures may result, paradoxically, in a generally stationary trend of total annual rainfall, in spite of the fact that daily rainfall increases [63]. In this context, it is relevant to learn how past warming has affected the changes of extreme precipitation [64,65].

To identify which monthly or seasonally rainfall could have led to a major increase in the annual storm erosivity in Southern Italy, we compared the rainfall intensity of the more recent period 1991–2014, with that one of the previous period 1961–1990. As shown in Figure 9, the most positive anomalies of rainfall intensity are observed mainly during September and December.

A positive core was extended from central Italy toward southern Italy in September, which might have triggered a high-velocity of soil loss, and heavy rain in the months in which the higher frequency of soil tillage was detected in the study area.

The observed increasing trend in the rain rates corresponds to the climate change prognosis that predicts reductions in average summer and autumn precipitation, combined with an increase in high-intensity rainy events for many parts of European agricultural lands [67,68]. Nearing et al. [69] point out that just because the sensitivity values for runoff and erosion were generally greater for rainfall changes as opposed to cover changes, this does not imply that future changes in rainfall will dominate over changes in land-use. Predictions by climate change experts suggest that the possibility of rainfall changes in the order of a few percent values in rainfall depths and intensities. In this context, with global and local climate change, an increasing number of ecosystems exist that are experiencing meteorological events of novel magnitudes, with timing and durations that are out of synchrony with the stress thresholds of organismic communities [70,71].

5. Conclusions

We have combined the RUSLE (Revised Universal Soil Loss Equation) and Thornes’ models to simulate gross soil erosion rates at the Met European Research Observatory of Monte Pino experimental site (Southern Italy), using plot-scale measurements and weather observations. The proposed method, resulting from a GIS-based upscaling procedure shows the fractal pattern that mirrors at a larger scale (sloped side surface) the soil loss data collected from plot-scale measurements. Model upscaling through scale-related parameters is thus a key to represent slope dynamics with plot-based, semi-empirical modeling. This approach proved effective to detect magnitude and timing of disproportional and nonlinear variations of soil erosion over a decadal time-scale in the Calore River Basin. Results from our modeling study indicate a clear increase in the peak values of gross erosion after 1995, with the highest peak value occurring in November 1997. Other peaks occurred in May 2001, and June 2003. Field observations revealed that it was mostly during long and moderately intense rainfall events occurring in winter that the sediment mobilized along the sloped side surfaces reached the main river, contributing to the sediment yield. These results thus suggest that, in the study area, the control exerted on soil erosion by spring and summer rainfall events of short duration and high intensity might be a less important part of the hydrological control on sediment yield. As revealed by some discrepancies between basin-wide estimates and measurements, hidden non-linear processes occur in down-and-up-scaling ranges [72], not easy to capture through parsimonious modeling. Presentation of these results is intended to provide an initial understanding for more detailed numerical modeling. For instance, the capacity to predict event soil losses could be improved by considering runoff explicitly [47,73]. This implies reporting discharge information or runoff data from the experimental plot, which were not available in this study.