Soil erosion of Hungary assessed by spatially explicit modelling

ABSTRACT

The conservation of soil resources is increasingly becoming a critical issue worldwide, with growing interest in carbon stocks and water storage within the soil. Hungary is no exception, and there has been a demand for a country level soil erosion map that incorporates digital information available from the latest surveys and digital mapping campaigns. The map presented in this paper in based on the extremely wet year of 2010, and thus provides users a 1:100,000 scale ‘worst case scenario’ of soil erosion risk in Hungary (see Main Map). Results from both the Universal Soil Loss Equation and the Pan-European Soil Erosion Risk Assessment models were combined in order to achieve a map that can be used by a wide range of professionals. Both models estimate soil erosion by water in tonnes per hectare per year.

KEYWORDS:

• Soil erosion
• PESERA
• USLE
• erosion model
• spatial risk assessment

1. Introduction

Soil erosion by water is a natural process driven, primarily by precipitation and topography. It is also influenced by other factors, such as soil characteristics, surface cover, temperature and land use practices. The European Commission recognized soil erosion as one of the main threats to soil resources in its Thematic Strategy on Soil Protection (2002).

There are a variety of approaches for mapping soil erosion susceptibility. The first widely used method for soil loss calculation was the Universal Soil Loss Equation (USLE) model created by Wischmeier and Smith (1978). Their work calls attention to the importance of local measurements and parametrization of the USLE model in order to reach proper values during the use of the model. The most important improvement to the model was the conversion of the dimensions from US to SI units (Foster, McCool, Renard, & Moldenhauer, 1981). Many scientists are using the USLE model for various purposes, even in situations, where in situ data are lacking. Baban and Wan Yusof (2001) used the USLE model to prepare an erosion risk map for Langkaey Island of Malaysia. Lahloi, Rhinane, Hilali, Lahssini, and Khalile (2015) used the USLE with remote sensing for potential erosion risk calculations in Morocco. USLE, and its revised version RUSLE, are used mainly for estimation of the yearly amount of soil loss in tonnes per hectare per year (t/ha/y) (Chen et al., 2010; Doğan et al., 2015; Fantappiè, Priori, & Costantini, 2014).

Input parameters of USLE and its developed versions MUSLE (Odongo, Onyando, Mutua, van Oel, & Becht, 2013; Sadeghi et al., 2007) and RUSLE (Conforti et al., 2015) have been thoroughly studied as summarized in Table 1.

Table 1. Studies on USLE input parameters.

Research	Country/Region	Publication
R factor	Sri Lanka	Joshua (1977)
	Belgium	Bolline et al. (1980)
	West Africa	Roose (1980)
	Southern Africa	Smithen and Schulze (1982)
	USA	Renard and Freimund (1994)
	Cape Verde	Mannaerts and Gabriels (2000)
	Brasil	da Silva (2004)
K factor	West Africa	Roose (1980)
	West Java	Ambar and Wiersum (1980)
	Brasil	Biscaia, Rufino, and Henklain (1981)
	Italy	Zanchi (1988)
	NA	Auerswald, Kainz, Angermüller, and Steindl (1996)
	China	Zhang, Shu, Xu, Yang, and Yu (2008)
	Greece	Gitas, Douros, Minakou, Silleos, and Karydas (2009)
LS	–	Desmet and Govers (1996), Formaggio, Gameiro, and Epiphanio (1998), Liu, Kiesel, Hörmann, and Fohrer (2011)
C	Spain	Folly, Bronsveld, and Clavaux (1996)
	Belgium	Gabriels, Ghekiere, Schiettecatte, and Rottiers (2003)
	Croatia	Basic, Kisic, Mesic, Nestroy, and Butorac (2004)
	Greece	Gitas et al. (2009)
	Spain	García-Ruiz (2010)
	Jordan	Alkharabsheh, Alexandridis, Bilas, Misopolinos, and Silleos, (2013)
	Greece	Alexandridis, Sotiropoulou, Bilas, Karapetsas, and Silleos (2014)
	Europe	Panagos et al. (2015a)
P	Texas, USA	Morgan and Nalepa (1982)
	–	Foster and Highfill (1983)
	Crete, Greece	Karydas, Sekuloska, and Silleos (2009)
	–	Chang, Hromadka, and Chang (2011)
	Europe	Panagos et al. (2015b)

As emphasized by Wischmeier and Smith (1978), in situ analyses of the modelled values should not normally be avoided. Many good examples for this practice exist: for example, Magliulo (2010) integrated the analysis of aerial photos and topographic maps with field observations. In the case of Hungary, however, there is a lack of a harmonized, country-wide erosion monitoring network that could support the appropriate validation of such model results (Kertész & Centeri, 2006). Such works, therefore, are limited to the use of relatively well-established models and the use of cross-validation between different modelling results.

The Pan-European Soil Erosion Risk Assessment (PESERA) model was developed for rill and inter-rill erosion estimates, primarily for applications at a regional scale, with a focus on Europe-wide input data availability (Kirkby, Irvine, Jones, & Govers, 2008). As opposed to the USLE, it is a process-based model, and has a significantly higher demand for input data, with a possible total of 128 data layers (Irvine & Kosmas, 2003). The model breaks up precipitation into overland flow, evapotranspiration and soil moisture. A generic model is used to calculate plant growth based on transpiration and land use information. The total amount of erosion is calculated based on soil erodibility, slope parameters and overland flow (Kirkby et al., 2008).

Hungary has a rich heritage of soil maps, with multiple mapping campaigns carried out over the twentieth century (Várallyay, 2002). However, the assessment of soil erosion was not the main purpose of these campaigns and therefore the first national level map of soil erosion was only published in the 1960s (Stefanovits, 1964).

With the advancement of computers and information systems, it later became possible to apply erosion models at a national scale. Centeri and Pataki (2000) applied the USLE at 1:100,000 scale for Hungary based on currently available spatial data layers.

Unfortunately, there is limited availability of continuous in situ monitoring data on soil erosion in Hungary (Kertész & Centeri, 2006). This limits the reliability of such models, but also raises the importance of cross-validation and the use of reliable, well-established models, such as the USLE and the PESERA models (Gobin, Govers, & Kirkby, 2006)

The aim of the work presented here was to develop an up-to-date soil erosion map for the whole area of Hungary, using a combination of the two methods discussed above, thus providing a common reference for future work in the field. There is a significant increase in the occurrence of extreme precipitation events in Europe in recent decades (Alexander et al., 2006; Klein Tank & Können, 2003; Zolina, Simmer, Belyaev, Kapala, & Gulev, 2009). Soil water erosion research should therefore focus on high precipitation years because they are when most of the severe soil erosion events will occur (van den Besselaar, Klein, Tank, & Buishand, 2012). Because it was a year with extreme precipitation rates (Bissolli, Friedrich, Rapp, & Ziese, 2011), 2010 was selected as the baseline year.

2. Materials and methods

2.1. USLE input data

In order to achieve compatibility between the two models, the same data sources were used as input data in both cases. The target cell size for the map grid was 100 m ×  100 m. Figure 1 presents the maps for the input factors used for the USLE model. The model area itself covers the whole of Hungary (approximately 93,000 km²).

Figure 1. Maps of USLE input factors.

Climatic information was derived primarily from the CARPATCLIM database (Szalai et al., 2013), with the addition of AGRI4CAST MARS data (https://ec.europa.eu/jrc/en/mars) where the former was not available. Where the data were not directly available from the dataset, it was calculated from the source data. For the USLE model, mean monthly rainfall data were used as a basis for calculating the R factor.

EU-DEM (Bashfeld & Keim, 2011) was used for the calculation of the terrain characteristics. In the USLE model, the LS factors describe the effect of topography on soil erosion, which combines slope length (L) and slope gradient (S). The LS factor was calculated using Moore’s method (Moore, Grayson, & Ladson, 1991) in SAGA GIS.

For land use related layers, the 100 m resolution CORINE Land Cover (CLC2006) dataset was used, with maize as a default crop cover for arable land.

Soil-related parameters were derived from the Digital Kreybig Soil Information System (DKSIS; Pásztor, Szabó, Bakacsi, Matus, & Laborczi, 2012) and Hungarian Soil Information and Monitoring System (SIMS; Várallyay, 2002) databases. Both USLE and PESERA models require different soil property maps to calculate the soil-related parameters. These maps were created using digital soil mapping (DSM) methods – these are spatially exhaustive, environmental auxiliary variables related to soil forming processes for the spatial inference of point data (McBratney, Mendonça-Santos, & Minasny, 2003). The most commonly used environmental covariates refer to topography, climate, lithology, land cover and land use. The soil-landscape relation can be modelled by geostatistical and data-mining methods. In the present study we used regression kriging, a method widely used in DSM. Regression kriging is a spatial prediction technique that combines the regression of the dependent variable on auxiliary variables with kriging of the regression residuals (Hengl, Heuvelink, & Rossiter, 2007).

In the USLE model, the K factor characterizes soil erodibility, which means the susceptibility of the soil to erosion by rainfall and runoff. Soil erodibility is primarily determined by particle size and organic matter content. A numeric estimation of K can be calculated by the following formula (Sharply & Williams, 1990), where SAN, SIL, CLA and OM are the sand, silt, clay and organic matter content of the soil and SN₁ = 1−SAN/100 (see also Figure 2).The following soil property maps were compiled to satisfy the demands of the USLE and PESERA models: particle size distribution maps (percentage sand, silt and clay), organic matter content and carbonate content according to the objectives of the DOSoReMI.hu initiative (Pásztor et al., 2015). The particle size maps were compiled using SIMS as a reference data source (see Laborczi, Szatmári, Takács, & Pásztor, 2015); the soil organic content and carbonate content maps were created based on the data of DKSIS.

Figure 2. Input data used for the calculation of the USLE K factor.

2.2. PESERA input data

Data sources for the PESERA model were the same as used for the USLE model input parameters, unless stated otherwise.

The PESERA model grid cell size had to be reduced to 120 m ×  120 m due to limitations in the model. The output was later resampled to match the USLE output grid.

The PESERA model required the following climate data for input: mean monthly rainfall, mean monthly rainfall per rain day, coefficient of variation of monthly rainfall per rain day, mean monthly temperature and monthly temperature range.

The topographic factor of PESERA focuses on the effect of local relief on soil erosion. The topographic factor is derived from the standard deviation of elevation within a 1.5 km distance around each cell.

For land use related layers, the 100 m resolution CORINE Land Cover (CLC2006) dataset was used, based on the PESERA manual (Irvine & Kosmas, 2003). In the case of arable land, maize was used as the crop cover.

In the PESERA model, the soil parameter is more complex and is composed of multiple factors. Soil erodibility is the most important factor, which characterizes the sensitivity of a soil to erosion. Soil erodibility is primarily related to soil texture, but it is also affected by organic matter and carbonate content. The erodible fraction (EF) was calculated by the following formula (Fryrear et al., 2000), where SAN, SIL, OM, CaCO₃ are the sand, silt, clay, organic matter and carbonate content and S_c is the ratio of sand to clay content.The mechanical impact of raindrops leads to soil surface crusting, which determines the lower limit of water storage capacity. The soil crust factor (SCF) is affected mostly by the clay content and organic matter content of the soil. This was calculated according to Fryrear et al. (2000):where CLA and OM are the clay content and the organic matter content.

The readily available soil water capacity provides the maximum storage capacity of the soil before runoff. PESERA derives this complex soil feature from three components: soil water available to plants, effective soil water capacity and scale depth. The applied Hungarian soil profile databases do not contain information on these parameters. To achieve the necessary input at the proper spatial resolution the original European PESERA data layers were downscaled. Virtual reference point sets were created by point sampling the low resolution maps. Ten points per 100 km² were randomized in the geographical space, representing virtual sampling locations. Conditional generalization was applied, prescribing a minimum spacing of 100 m (equal to the cell size applied in spatial modelling) between the generated points. The values of the dependent (predicted) and independent (predictor) variables were identified at the randomized locations and their records were used in the further spatial inference. Since all three parameters were quantitative, regression kriging was also applied in the disaggregation process.

3. Results

3.1. Running the models

Both the USLE and the PESERA models have been applied for the year 2010. Due to the extreme precipitation rates of that year, the output layers of both models represent a ‘worst case scenario’ from a climatic point of view for soil erosion potential in Hungary. The extreme precipitation represents a 20 year return frequency. The models are limited by a lack of information for arable land. We can only prepare scenarios for a given crop, crop rotation or land cover for arable areas as the exact crop rotation information is not available at a national scale. In this present case, the C factor of maize was applied, which represents a ‘worst case scenario’.

Raw, unclassified output results from the USLE and PESERA models are presented in Figures 3 and 4, respectively. Since both methods provided the output as t/ha/y, soil loss comparison and combination of the two maps is possible.

Figure 3. Output results from the USLE model.

Figure 4. Output results from the PESERA model.

3.2. Harmonization of model results

In order to provide one unified map of results, a simple ensemble model was applied. The mean value of the two output maps was calculated for each grid cell. Detailed analysis of the differences between results is beyond the scope of the present work and will be published as a separate study.

Figure 5 presents the areas affected by soil erosion as a percentage of Hungary’s total land area. About 74% of the area is affected only weakly, 18% is moderately affected, and only about 8% is affected by strong erosion.

Figure 5. Spatial extent of erosion categories in Hungary (as % of total area).

The resulting map – the main result of the present study (see Main Map) – provides an overview not only of the rate of erosion in Hungary, but also its spatial distribution. Looking at the map it is clear how the values presented in Figure 5 relate to Hungary’s geography. As expected, high erosion rates mostly occur in mountainous and hilly regions, such as Northern Hungary, or Southern Transdanubia. The latter shows the highest erosion rates in some areas. This is due to the fact that the mountains in the north have greater forest coverage, while in the south many hills are cultivated. The Great Plain in general shows low erosion rates, with a few regions where mostly sandy soils dominate. The large area covered by the plains explains the high percentage of area affected only by low erosion rates, while it is clear that in the hills and mountains erosion rates are significantly higher.

4. Conclusions

The map presented here is the first attempt in Hungary to provide a model-based soil erosion map at the national level, combining two previously applied methods. The use of a harmonized input data set provided a base for comparison and synthesis of the two methods.

The results show that about 26% of the country’s total area is prone to moderate to high erosion rates. The map presented here supports previous studies, presenting a distribution of erosion risk that is similar to earlier works in both spatial patterns and magnitude. The new map is more detailed and the use of the best available data sets also makes it more reliable. Nevertheless, it should be noted that in order to provide a more dependable output, field validation of the model results would be required. To be able to do that, a country-wide network of erosion monitoring sites is needed as this does not presently exist in Hungary.

Future research will focus on a more detailed evaluation of combining the two models, with special focus on field validation of the modelled results, and on local discrepancies.

Software

SAGA GIS 2.2.0 (Böhner, McCloy, & Strobl, 2006) was used for the preparation of data layers. Esri ArcGIS 9.0 and ArcInfo Workstation 9.0 were used for running the PESERA model, including pre- and post-processing. Spatial analysis and map layout compilation were carried out in Esri ArcGIS 10.3.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

Our work has been supported by the Hungarian National Scientific Research 285 Foundation (OTKA) [grant number K105167) and the Hungarian grant Research Centre of Excellence-11476-3/2016/FEKUT.