Towards the Assessment of Soil-Erosion-Related C-Factor on European Scale Using Google Earth Engine and Sentinel-2 Images

Abstract

Soil erosion is a constant environmental threat for the entirety of Europe. Numerous studies have been published during the last years concerning assessing soil erosion utilising Remote Sensing (RS) and Geographic Information Systems (GIS). Such studies commonly employ empirical erosion models to estimate soil loss on various spatial scales. In this context, empirical models have been highlighted as major approaches to estimate soil loss on various spatial scales. Most of these models analyse environmental factors representing soil-erosion-influencing conditions such as the climate, topography, soil regime, and surface vegetation coverage. In this study, the Google Earth Engine (GEE) cloud computing platform and Sentinel-2 satellite imagery data have been combined to assess the vegetation-coverage-related factor known as cover management factor (C-factor) at a high spatial resolution (10 m) considering a total of 38 European countries. Based on the employment of the RS derivative of the Normalised Difference Vegetation Index (NDVI) for January and December 2019, a C-factor map was generated due to mean annual estimation. National values were then calculated in terms of different types of agricultural land cover classes. Furthermore, the European C-factor (C_EUROPE) values concerning the island of Crete (Greece) were compared with relevant values estimated for the island (C_CRETE) based on Sentinel-2 images being individually selected at a monthly time-step of 2019 to generate a series of 12 maps for the C-factor in Crete. Our results yielded identical C-factor values for the different approaches. The outcomes denote GEE’s high analytic and processing abilities to analyse massive quantities of data that can provide efficient digital products for soil-erosion-related studies.

1. Introduction

Soil erosion is a physical process that constitutes a major environmental threat in the European continent and worldwide [1]. xtensive research initiatives have been applied to support the study of soil erosion monitoring, assessment and mitigation in recent years [2,3]. Geospatial technologies such as Remote Sensing (RS) and Geographic Information Systems (GIS) have played a crucial role in this task, offering the opportunity to researchers to explore extensive regions on various time and spatial scales. Multiple modelling approaches have been developed and utilised for soil erosion assessment. Empirical models can be characterised as the most widely applied and accepted worldwide [4]. The Universal Soil Loss Equation (USLE) and its revised version (Revised Universal Soil Loss Equation—RUSLE) constitute the main representatives of empirical models [5].

The aforementioned empirical models are commonly based on integrating natural and human-induced environmental factors describing the main components of the erosion process. One of them is the cover management factor (C-factor), which can be expressed as a soil loss relation of a given plot of land covered with specified vegetation to a bare seedbed-prepared plot ploughed up and down along the slope gradient [6,7,8]. Being affected by rainfall conditions and human interventions in a region, C-factor is negatively related to soil loss rates [9]. To estimate the C-factor, different approaches have been employed, such as applying uniform values from literature to land use/cover classes [10] or generating vegetation indices such as the Normalised Difference Vegetation Index (NDVI) [11]. Usually, NDVI is directly correlated to the C-factor by a linear or exponential regression [12]. However, previous studies have observed that NDVI-based estimation of the C-factor is usually sensitive to various biophysical variables, such as soil regime, topographic features, and vegetation phenology [13]. Thus, discrepancies are usually identified between the NDVI-based and literature-based factor values [14].

Google Earth Engine (GEE) is a cloud computing platform that specialises in processing satellite imagery and provides tremendous image data support. This high-performance computing infrastructure enables researchers to easily and quickly access more than thirty years of freely distributed public archives for developing RS-based applications [15]. GEE is open, easy to develop algorithms, and can batch process image data quickly, reducing the cost and complexity of geospatial data analysis compared with traditional image processing tools. Through GEE, algorithms and results are freely shared and validated. The function of mass data storage, which can quickly invoke and batch process huge amounts of data, makes GEE an ideal tool to assess and map phenomena on national, continental or even global scales.

The main objective of this study is to explore the capabilities of GEE to effectively estimate C-factor related to land cover on a European scale by using NDVI data products from Sentinel-2 satellite images. C-factor values were estimated and mapped for 38 European countries, as well as for their different agricultural lands represented by different land cover classes. In this context, we mainly focused on agricultural and pasture land classes since (a) the agricultural areas are more prone to soil erosion, (b) the NDVI approach is more efficient on these classes and (c) the soil loss phenomena have an economic impact on these land use classes, thus providing valuable information. Additionally, in order to estimate the accuracy of the proposed methodological framework, the results of the European-scale assessment were compared with the results derived both from relevant research works [10,16] as well as from the regional-scale assessment of C-factor values as estimated in the island of Crete, Greece.

2. Materials and Methods

This study was conducted for 38 countries of Europe, covering the continent’s main part from Portugal to the west, the European part of Russia to the east, Norway to the north and Greece to the south. Thus, a total amount of 36,936 Sentinel-2 satellite images, covering the period between 1 January 2019 and 31 December 2019, were analysed to estimate the C-factor across Europe. The Sentinel-2 Multi-Spectral Instrument (S2-MSI) mission was launched by the European Space Agency (ESA) in 2015. Two satellites (2A and 2B) were placed in orbit, which combined provided an approximately five-day revisit time. Sentinel 2-MSI has spatial resolutions ranging from 10 to 60 m and thirteen spectral bands that cover the visible (coastal, blue, green and red), near-infrared (B8 and NIR), red-edge (RedEdge1, RedEdge2, and RedEdge3), shortwave-infrared (SWIR), and water and cirrus regions [17]. The overall methodology applied in the GEE environment is presented in Figure 1.

A significant number of Sentinel-2 images (10 m spatial resolution) were employed in the GEE platform to create an overall NDVI-based C-factor mosaic for Europe. Crucial parameters for accomplishing the specific task were determining a universal cloud coverage threshold of <20% and the atmospheric correction for all the images to be used. After selecting, semi-automatically, the appropriate images, the overall mosaic was developed. At the next step, the NDVI was calculated as follows:

$$NDVI=\frac{\left(NIR-RED\right)}{\left(NIR+RED\right)}$$

(1)

In cases of images overlapping, the mean NDVI value was considered. Then, the C-factor was estimated according to the following equation [17]. According to [15], the NDVI exponential function gives more accurate results than the linear regression NDVI approach that describes +0. The linear relationship between the NDVI and the C-factor is as follows:

$$C=exp\left[-a*\frac{NDVI}{\left(b-NDVI\right)}\right]$$

(2)

where a and b are fitting parameters with values 2 and 1, respectively. This equation has been used in numerous research efforts to generate C-factor [18,19,20], some of them applied in the European continent [18], making the choice of the specific values ideal for our study. The estimated C-values were then cartographically visualised as a map to show the spatial distribution of the factor across Europe (Figure 2). The vegetation coverage and, consequently, the corresponding NDVI values are inversely proportional to soil erosion phenomena, which means that areas with high NDVI values are less prone to soil erosion (low C-values). Furthermore, topographical variations also affect the NDVI-derived C-values by affecting both the (a) spectrum reflectance property of the surface and (b) the greenness of the vegetation [19].

In terms of further processing, a random dense point sampling (5000 points for each county) was carried out over the European C-factor product to estimate the mean C-value agricultural land cover classes defined by the “CORINE” 2018 database (

Table 1. Description of agricultural land cover classes, according to “CORINE” database.

CORINE Land Cover Code	Description
211	Non-irrigated arable land
212	Permanently irrigated land
221	Vineyards
222	Fruit trees and berry plantations
223	Olive groves
231	Pastures
241	Annual crops associated with permanent crops
242	Complex cultivation patterns
243	Agriculture areas with significant areas of natural vegetation

and

Table 2. Mean C-factor values for the European countries and their various agricultural land cover classes.

Country	Overall	CORINE Land Cover Code
		211	212	221	222	223	231	242	243
Albania	0.145	0.197	0.125	0.118	0.098	0.083	0.22	0.16	0.157
Austria	0.107	0.184	-	0.172	-	-	0.076	0.053	0.052
Belgium	0.064	0.136	-	-	0.044	-	0.034	0.059	0.046
Bosnia & Herzogovina	0.071	0.131	0.061		0.044	-	0.083	0.058	0.051
Bulgaria	0.185	0.301		0.212	0.148	-	0.175	0.126	0.119
Croatia	0.116	0.157	0.175	0.158	-	0.148	0.085	0.091	0.059
Cyprus	0.389	0.495	0.44	0.376	0.292	0.421	0.288	0.411	0.369
Czech Republic	0.128	0.171	-	0.211	0.123	-	0.091	0.083	0.086
Denmark	0.07	0.089	-	-	0.088	-	0.055	0.057	0.06
Estonia	0.146	0.19	-	-	0.174	-	0.135	0.131	0.102
Finland	0.274	0.292	-	-		-	0.01	0.191	0.26
France	0.13	0.165	-	0.202	0.094	-	0.06	0.083	0.068
Germany	0.073	0.123	-	0.094	0.052	-	0.056	0.04	0.071
Greece	0.204	0.299	0.212	0.184	0.101	0.156	0.232	0.216	0.156
Hungary	0.172	0.212	-	0.196	0.087		0.12	0.119	0.086
Iceland	0.154	0.04	-	-	-	-	0.208	0.215	-
Ireland	0.054	0.11	-	-	-	-	0.019	0.06	0.028
Italy	0.179	0.237	0.224	0.234	0.111	0.168	0.151	0.176	0.104
Latvia	0.141	0.151	-	-	0.228	-	0.102	0.131	0.09
Lithuania	0.108	0.152	-	-	0.085	-	0.091	0.127	0.084
Luxembourg	0.062	0.089	-	0.087	-	-	0.031	0.057	0.045
Malta	0.373	0.363	-	0.338	-	-	0.377	0.411	0.381
Montenegro	0.093	0.108	-	-	0.067	0.021	0.18	0.109	0.087
Netherlands	0.075	0.169	-	-	0.072	-	0.032	0.071	0.033
North Macedonia	0.255	0.35	0.321	0.324	0.16	-	0.276	0.278	0.172
Norway	0.166	0.238	-	-	-	-	0.054	0.169	0.204
Poland	0.103	0.162	-	-	0.062	-	0.075	0.128	0.09
Portugal	0.221	0.302	0.154	0.254	0.232	0.245	0.262	0.219	0.169
Romania	0.154	0.243	0.137	0.177	0.073	-	0.124	0.135	0.084
Serbia	0.165	0.225	-	0.219	0.137	-	0.156	0.151	0.098
Slovakia	0.132	0.2	-	0.094	0.14	-	0.079	0.19	0.086
Slovenia	0.067	0.114	-	0.045	0.118	0.062	0.043	0.053	0.037
Spain	0.4	0.463	0.372	0.353	0.384	0.4	0.215	0.267	0.317
Sweden	0.157	0.164	-	-	-	-	0.158	0.182	0.123
Switzerland	0.101	0.059	-	0.098	0.174	-	0.119	0.047	0.111
United Kingdom	0.049	0.124	-	-	0.018	-	0.03	0.046	0.027
Europe	-	0.2	0.238	0.206	0.126	0.189	0.118	0.132	0.114

) on a national level. The “CORINE” database is the de facto standard option for land use/cover (LULC) monitoring at the European level. It consists of an inter-annual LULC dataset for Europe, produced by national agencies and coordinated by the European Environment Agency (EEA) [21]. This paper focuses mainly on arable land classes because: (a) the NDVI approach can give more efficient results in these areas and (b) the soil loss has a negative effect on agricultural areas, and this may provide valuable info on a European scale.

3. Results and Discussion

The overall C-factor results denote that the higher mean values for 2019 are observed in the island countries of Europe and specifically Cyprus and Malta. Furthermore, except for the case of Finland, where quite high values are indicated (maybe this is a drawback due to the presence of snow cover), countries of southern Europe such as Spain, Portugal, North Macedonia and Greece seem to also be prone to soil erosion due to high corresponding C-factor values. On the other hand, countries of northern Europe such as Iceland, the UK and Belgium have generally low values ranging around 0.05.

C-factor values were also estimated concerning the “CORINE” agricultural classes for each country. It is worth mentioning that particular focus was given on arable land classes (codes of “211” and “212”) where the soil erosion process has higher economic effects. In this vein, high C-values were indicated for the “non-irrigated arable land” class (code of “211”) in countries of southern Europe such as Greece, Cyprus and Malta rather than central Europe (e.g., Switzerland, Poland), where values were found to range around 0.04. Regarding the rest of the agricultural classes, a similar (to arable land) response was revealed for “vineyards” and “complex cultivation patterns” (codes of “221” and “242”, respectively) where, again, southern Europe presents greater values than central and northern Europe. The overall results concerning the estimation of C-factor in Europe were compared to the results of [10], where both land cover and management practice parameters, as estimated from literature values, were incorporated in the overall analysis. Although country-level differences may be observed between the two methodologies, the overall mean values in terms of Corine land cover classes are quite similar. Thus, for our study, the mean C-value compared to the results of [10] is: (a) 0.20 vs. 0.27 for class 211; (b) 0.238 vs. 0.298 for class 212; (c) 0.2069 vs. 0.298 for class 221; (d) 0.126 vs. 0.20 for class 222; (e) 0.1892 vs. 0.208 for class 223; (f) 0.118 vs. 0.111 for class 231; (g) 0.15 vs. 0.139 for class 242; and (h) 0.114 vs. 0.119 for class 243. Although the estimation of NDVI-based C-factor has specific limitations, such as the lack of incorporation of management practices in the overall methodology, it manages to give emphasis on the temporal variation of NDVI sensitivity and contributes to the calculation of C-factor in extended areas.

Comparison to Regional Assessment

In order to assess the overall credibility of the applied methodology for European-scale assessment, the island of Crete in Greece—the largest Greek island and fifth-largest Mediterranean island—was used as a pilot area for larger-scale (regional) assessment. In particular, 12 Sentinel-2 satellite images were “handpicked” and analysed in the GEE environment for NDVI (Equation (1)), and thus C-factor estimation (Equation (2)) at the monthly time-step of 2019 (Figure 3) was estimated. The main criterion of imagery selection was the low cloud cover percentage of the images in the pilot area. In this context, Crete was used as a “controlled” pilot area, where a single image for every single month of 2019 was selected, avoiding the random incorporation of numerous images as in the European example.

The C-factor results derived from the Cretan assessment (C_Crete) were compared with those from the European assessment (C_Europe). Evaluation of regression fit was carried out using Pearson’s correlation coefficient to investigate the correlation between C_Europe and C_Crete values. This coefficient is a measure of linear dependence between two variables, x_i and y_i, with values ranging between [−1 +1] (a value of 1 if the two sample variables are identical, and a value of -1 if they are completely uncorrelated) [22]. The Pearson correlation coefficient between C_Europe and C_Crete (annual mean value) was found to be around 0.7, indicating a moderate-to-high correlation between the two assessments of different scales. Regarding their response to the examined land cover classes, the C_Europe and C_Crete values were generally consistent, with relatively minimal differences between them. Specifically, the two-factor assessments gave almost identical results for the classes coded by “212”, “223”, “242”, and “243” (Figure 4). Furthermore, their differences concerning classes of “211”, “221”, “222”, and “231” could be considered to be within the standard deviation error, related to the fact that these classes cover smaller areas in the island. All of the above revealed that, in the NDVI-based approach, the semi-automatic development of a continental C-factor mosaic offered more or less similar results to a regional factor product that was developed using only selected Sentinel-2 images.

In addition, the results of both approaches seem to agree with the Greek C-factor results estimated by a literature-based approach in the work of Panagos et al. [10] and Karydas et al. [16]. According to the specific works, the mean C-factor value (0.28) for arable lands in Greece and Crete, respectively, was very close to both C_Europe and C_Crete values (0.26 and 0.24, respectively). Although the spatial resolution (100 m) in the analysis of [8] is much higher than this one of Sentinel-2 imagery data (10 m) used in this study, the comparison can provide answers to the common debate between the efficiency of either NDVI-based or literature-based approach for C-factor assessment. In general, the NDVI-based approach overcomes the main disadvantage of the literature-based approach, which usually incorporates static C-factor values ignoring the chronological and spatial dynamics of land use components [23,24]. Therefore, approaches that use Earth Observation (EO) data in the cloud computing environment appear as an alternative to minimise this limitation and provide freely available data of an adequate spatio-temporal resolution.

4. Conclusions

Empirical modelling is the most widely used approach for soil erosion assessment and soil and water conservation planning. One of the most crucial parameters in the majority of empirical models is the soil-erosion-related C-factor. This study presents a sophisticated way to assess and map the C-factor within the GEE cloud computing environment using Sentinel-2 imagery data.

From a general perspective, the overall methodology forms a sufficient framework for assessments on coarse scales such as those of the European continent but based on detailed, high-quality data (with a spatial resolution of 10 m). However, it has to be pointed out that the NDVI approach has specific limitations, such as that it focuses only on land and crop cover and not on soil management. The GEE’s ability to analyse a massive amount of data (more than 35,000 satellite images) and provide efficient digital products with minimal computer time processing consumption is highlighted.

Furthermore, the potential of the synergistic use of Sentinel-2 imagery data and the GEE environment is pointed out to provide fast and efficient results for soil erosion monitoring.

Future research work will focus on further in situ soil loss data collection for more accurate validation of this study’s findings and the production of additional erosion-influencing factors within the GEE environment to be incorporated in empirical modelling approaches.