The Use of Spectral Indices to Recognize Waterlogged Agricultural Land in South Moravia, Czech Republic

Abstract

The agricultural landscape of the Czech Republic is facing climate change, and drought is among the most severe stress factors. Thousands of small ponds and naturally wet areas have been drained and transformed into agricultural parcels. Their restoration could increase the landscape’s resilience to climate change. Therefore, we describe the possibility of using hyperspectral aerial surveying for the identification of waterlogged areas in the agricultural landscape based on the example of one of the warmest and driest regions of the Czech Republic—the South Moravian region, an area where water retention in the landscape is highly relevant. Within our study, a total of 33 spectral indices related to the waterlogging of soil selected from previous studies were evaluated. The maximum entropy model (MAXENT) was used in the analysis of these indices. The analysis, which was carried out in several locations during different periods of the year (spring and autumn), shows the varying applicability of individual groups of indices. Regardless of the season, chlorophyll-based indices (MCARI—31.8, CARI—26.3, TCARI2—24.3 average percentage contribution) made the most significant contribution to the creation of probability maps of the occurrence of waterlogged areas. However, more accurate results could be achieved in the spring period by using the NVI index (40.5 average percentage contribution). The results show that remote sensing could be used for the identification of waterlogged sites, especially for initial identification, which should then be confirmed by field survey. Furthermore, the research points out the role of the LAI and chlorophyll content. According to the NVI, low LAI contributes the most to the probability of occurrence in the spring season, while chlorophyll-based indices prove to be the best, contributing high values, which is rather contradictory but could be resolved only by subsequent field research.

1. Introduction

The issue of drought is much debated, particularly in recent years, in relation to climate change. Although total precipitation is decreasing in the Czech Republic, the timing of precipitation plays a large role, causing so-called agriculture drought, when crops do not have a sufficient amount of water in a given period of time [1]. Therefore, identifying waterlogged areas is important for appropriate measures to retain water and make it available during dry periods. For this reason, discussion also focuses on storm rainfall with consequential erosion processes and flooding [2,3,4].

In the current agricultural landscape of the Czech Republic, we can encounter many types of water and waterlogged areas that can be used as elements of landscape water management. In the past, small bodies of water were created in the agricultural landscape, mainly in the form of ponds, but later most often as multipurpose reservoirs, e.g., with irrigation or flood control functions. Today, these are also found in the form of (semi)dry polders, which are a specific case of (semi)dry reservoirs. The ponds were a typical component of the Czech landscape, built from the early Middle Ages, but they have decreased in number over the last few centuries from about 75,000 to 25,000 and have mainly been converted into agricultural land. In the period from the mid-19th century alone, over 3400 defunct ponds larger than 0.5 ha have been mapped [5,6,7]. We can still find traces of them in the local landscape today, for example, in the form of remnants of dikes or waterlogged sites that make farming difficult (Figure 1 and Figure 2). These can be areas where, under certain circumstances, we can consider restoring these bodies of water.

Waterlogged localities are also often related to disturbed drainage systems, which in the Czech countryside cover roughly 25% of the total agricultural land and were usually built more than thirty years ago. It is not uncommon to find drainage systems which are more than 100 years old. Within these, the most common defects are caused by silting, root penetration, mechanical damage, corrosion etc. [6,7]. On waterlogged sites with defective drainage, we must consider resolving these problems, while taking into account water retention in the landscape. This may entail not only restoring these facilities, but the possibility of using certain locations as wetlands or converting these facilities into two-way controlled systems capable of either draining water or retaining it in the soil. Or, if necessary, supplying it for irrigation or replenishing groundwater.

Waterlogged sites also include areas of agricultural land with signs of waterlogging (hygrophilous vegetation, stagnant surface water, etc.), with specific pedological and hydropedological properties, making it impossible to use these areas at least for part of the season. These also include periodically flooded ephemeral pools on arable land, which serve as a refuge for many wetland organisms in an otherwise homogeneous agricultural landscape and which host important communities of endangered and rare protected species of plants [8] and animals [9]. These organisms are bound to open habitats, which are usually in early stages of succession due to regular disturbance [10].

Within landscape planning, it is necessary to work with all these elements in the agricultural landscape and systematically implement a complex of measures that will be effective in terms of drought, soil conservation, and increased biodiversity. However, first of all, it is necessary to identify these areas in the landscape. Knowledge of the spatial and temporal occurrence of waterlogged areas is essential for their protection and development. For the mapping of these specific retention elements in the landscape, it is advisable for many reasons to use GIS data and tools or remote sensing. The development and availability of these geospatial technologies have recently led, in particular, to much better spatial, temporal, and spectral resolution of the image data used [11,12]. Thanks to this, it is currently possible to use a whole range of new tools when processing these data, e.g., in the form of remote sensing and spectral indices.

Spectral indices are generally used to increase the sensitivity of identification of the studied phenomenon to separate spectral bands or reflectance in individual wavelengths [13]. Vegetation indices are a subgroup of spectral indices which mainly show vegetation characteristics but also the soil properties which are often linked to them. In principle, we can distinguish between broadband indices (calculated from so-called multispectral data) which capture a wider range of the spectral band and narrow-band indices with precisely defined wavelengths calculated from so-called hyperspectral data [14,15].

Broadband vegetation indices particularly show the vigour and health of green vegetation. They are sensitive to the combined effect of leaf cover and grouping and chlorophyll content. In particular, the difference between the strong reflectance of vegetation in the near-infrared part and the low reflectance in the spectral part of red light is used. Compared to wide-band indices, narrow-band indices are more sensitive to minor changes in vegetation parameters, especially where vegetation density plays a role when narrow-spectral images tend to be oversaturated [16].

The potential of narrow-spectral indices for the identification of some vegetation parameters is known, but their targeted use for the analysis of waterlogged areas on agricultural land has not yet been tested more widely in the conditions of Central Europe.

Several approaches can be used to evaluate the suitability of spectral indices for identifying waterlogged areas. Statistical methods based on logistic regression, GLM [17], or even more complex nonparametric methods in the form of the random forest approach are used [18]. Methods based on machine learning, in the form of neural networks [19], or, e.g., genetic algorithms (GARP [20]) or heuristic models (BIOCLIM [21]) are also increasingly used. One of the models that have been successfully applied in various areas of interest is the maximum entropy model (MAXENT), used especially in the spatial distribution of species [22]. In principle, however, it is suitable for evaluating the occurrence of any phenomenon where environmental conditions become variables that cover the entire extent of the studied area. It makes it possible to create maps of probable occurrence, and also enables the automatic evaluation of the influence of variables in the entire range of their values. It works only with presence data, when it is not assumed that these data will be complete. In the case of waterlogged surfaces, this input set is somewhat uncertain, but this does not represent a major problem for the model.

Our study aimed to evaluate the most favourable spectral indices for identifying waterlogged areas on agricultural land which would enable us to quickly locate suitable areas and propose appropriate measures for water retention and increasing biodiversity in the agricultural landscape.

2. Materials and Methods

2.1. Study Area

The area of interest is located in South Moravia, SE of Brno, the extent of which corresponds to aerial hyperspectral images taken in 2016, 2017, and 2020 by the Flying Laboratory of Imaging Systems (FLIS) operated by the Global Change Research Institute.

From the geomorphological perspective, the area comprises mainly plains and hills. Elevation ranges approximately between 200 and 400 m above sea level. In the western part, the land rises from SW to NE, in the eastern part, it rises and then falls again in the NW–SE axis. The area is within a very warm-to-warm climate zone. At lower elevations, there is a typical seasonal lack of precipitation. At lower elevations, there are soils of chernozem groups, while at higher elevations, there are leptosols, luvisols, and cambisols. Around water courses, there are also fluvisols. Hydrologically, there are several subbasins which form the source area, the upper and middle parts of local water courses that flow into the Black Sea basin.

The area of interest is comprised of heavily used agricultural lands with a high proportion of fertile land. Wheat, barley, corn, and winter rape are the main crops grown, while the use of perennial fodder is low due to a reduction in livestock production. The area also includes the wine region of South Moravia. There are smaller industrial and food businesses in the area, and residents often commute to the nearby regional city of Brno for job opportunities.

A total of four locations were considered, three within one time period and one in four time periods (Figure 3). The area under consideration is an intensively farmed area with a predominance of chernozem soils in the warmest and driest climatic region of the Czech Republic. At the time of the aerial imaging campaigns, data were available on the relative saturation of the soil and the intensity of drought, as well as monthly precipitation totals in the monitored years. These data are shown in

Table 1. Soil wetness of individual locations about one week before aerial imaging (INTERSUCHO. Available online: https://www.intersucho.cz/ (accessed on 12 December 2022)).

Location	Campaign	Shooting Time	Relative Water Saturation of Soil *	Drought Intensity **
1	30 April 2016	10:43–11:53	70 ± 20	0
1	27 August 2016	10:36–12:34	55 ± 5	0
1	19 May 2017	14:13–15:09	65 ± 25	0
1	30 August 2017	11:20–13:34	35 ± 5	0
2	21 August 2020	10:14–11:31	55 ± 10	0
3	21 August 2020	10:14–11:31	50 ± 10	0
4	21 August 2020	10:14–11:31	45 ± 5	0

* Percentage of water saturation of the soil layer (0–40 cm). ** Deviation of soil moisture (expressed as the degree of drought) from the usual state from 1961 to 2010 in the 0–40 cm soil layer (range, 0–6, where 0—best, 6—worst).

and Figure 4.

2.2. Waterlogged Areas—Specification and Data Source

The basic underlying data layer for the vectorization of waterlogged areas and subsequent analysis of spectral indices was a seamless orthophoto map of the Czech Republic from 2000 to 2020. This data layer is publicly available on the ČÚZK geoportal (Czech Land Survey and Cadastral Office) as a WMS or viewing service. A total of nine time periods were used. Specifically, these were archival orthophoto maps from 2000, 2003, 2006, 2009, 2012, 2014, 2016, 2018, and the current orthophoto from 2020. The target waterlogged areas were first visually identified within the scope of the area of interest and on individual orthophoto maps (periods) and then continuously manually vectorized in the ArcGIS software environment.

Figure 5 shows typical waterlogged areas on agricultural land at different stages of development. Figure 5a shows a waterlogged area (red polygon) on part of a soil block in the stage of developed field wetland with hygrophilous vegetation; on the surrounding soil block, there was a cereal crop. In this image, a waterlogged area is also found in the flooded area of a defunct historical pond (see the blue border). In Figure 5b, the waterlogged site is in the initial stage of field wetland development with residual stagnant water on the surface without wetland vegetation. On the surrounding area of the soil block, there is oilseed rape just before flowering. On orthophoto images, it is sometimes possible to also observe, besides the waterlogging itself, its cause, which is often non-functional subsurface drainage (Figure 5d), which shows the same location as in Figure 5c, but in a different year. In the central part, the influence of subsurface drainage in the form of a series of collection and drainage channels is evident in vegetation growth. At the same time, when comparing the local situation, the development of the waterlogged area over the years as a result of non-functional drainage can be seen.

2.3. Hyperspectral Data

The spectral indices were computed from hyperspectral airborne images acquired with CASI and SASI spectroradiometers (Itres, Ltd.) on the Flying Laboratory of Imaging Systems (FLIS) operated by CzechGlobe on 30 April 2016, 27 August 2016, 19 May 2017, and 30 August 2017 (Figure 3, Locality 1), and then on 21 August 2020, in areas 2, 3, and 4 (see Figure 3).

Firstly, all images were radiometrically corrected, georeferenced, corrected for atmospheric influence, and combined into one data cube (185 spectral bands in steps of 371 and 2442 nm, with a bandwidth of 7 nm in the VNIR range and 15 nm in the SWIR range) [23]. Combination into one data cube required (Figure 6) downsampling of the spatial resolution of CASI images to match the nominal spatial resolution of SASI images (5.0 m). Apart from the actual spectral indices, the images were classified in ENVI utilizing the maximum likelihood supervised classification, differentiating between vegetation cover, bare soil, and other forms of land use.

2.4. Spectral Indices

Based on professional literature [13,24,25], 33 indices were selected that cover the broadest possible range of vegetation and soil manifestations of waterlogged areas (

Table 2. Spectral indices considered for the identification of waterlogged areas on agricultural land. Chl—chlorophyll, LAI—leaf area index, SM—soil moisture, WC—water content, DMC—dry matter content.

Spectral Index	Formula	Scale	Sensitivity	Reference(s)
CARI	$\left(R_{700}-R_{670}\right)/0.2\left(R_{700}+R_{670}\right)$	Leaf	Chl	[25]
DVI	$R_{800}-R_{680}$	Leaf, soil	WC	[26]
EVI	$2.5(R_{800}-R_{670})/\left(R_{800}-6R_{670}-7.5R_{475}+1\right)$	Canopy	LAI	[27]
GNDVI	$\left(R_{750}-R_{550}\right)/\left(R_{750}+R_{550}\right)$	Canopy	Chl a	[28]
MCARI	$\left[\left(R_{700}-R_{670}\right)-0.2\left(R_{700}-R_{550}\right)\left(R_{700}/ R_{670}\right)\right]$	Canopy	Chl, LAI	[29]
MCARI/ OSAVI	MCARI/OSAVI	Canopy	Chl, LAI	[30]
mNDVI705	$\left(R_{750}-R_{705}\right)/\left(R_{750}+R_{705}+2 R_{445}\right)$	Leaf	Chl, LAI, WC	[31]
MSAVI	$0.5(2R_{800}+1-{[{\left(2R_{800}+1\right)}^2-8\left((R_{800}-R_{670}\right)]}^{0.5}$	Canopy	Chl	[32]
MTVI1	$1.2\left[1.2\left(R_{800}-R_{550}\right)-2.5\left(R_{670}-R_{550}\right)\right]$	Canopy	LAI	[30]
MTVI2	$\frac{1.5\left[1.2\left(R_{800}-R_{550}\right)-2.5\left(R_{670}-R_{550}\right)\right]}{[{\left(2R_{800}+1{)}^2-\left(6R_{800}-5\sqrt{R_{670}}\right)-0.5\right]}^{0.5}}$	Canopy	LAI	[28]
NDCI	$\left(R_{762}-R_{527}\right)/\left(R_{762}+R_{527}\right)$	Lakes	Chl a	[33]
NDDI	$(\mathit{NDVI}-\mathit{NDMI})/(\mathit{NDVI}+\mathit{NDMI})$	Soil	SM	[34]
NDMI	$\left(R_{860}-R_{1240}\right)/\left(R_{860}+R_{1240}\right)$	Soil	DMC	[35]
NDRE	$\left(R_{740}-R_{705}\right)/\left(R_{740}+R_{705}\right)$	Canopy	Chl, LAI	[36]
NDVI	$\left(R_{800}-R_{680}\right)/\left(R_{800}+R_{680}\right)$	Canopy	Chl, LAI	[37]
NDVI705	$\left(R_{750}-R_{705}\right)/\left(R_{750}+R_{705}\right)$	Leaf	Chl a	[38]
NDWI	$\left(R_{560}-R_{860}\right)/\left(R_{560}+R_{860}\right)$	Water bodies, canopy, soil	WC	[39]
NMDI	$\left(R_{860}-R_{1240}\right)/\left(R_{860}+R_{1240}\right)$	Soil, canopy	SM, WC	[40]
NVI	$(R_{777}-R_{747})/ R_{673}$	Leaf	LAI	[41]
OSAVI	$[\left(1+0.16\right)\left(R_{800}-R_{670}\right)/\left(R_{800}+R_{670}+0.16\right)$	Canopy, soil	Chl, SM	[42]
PRI	$\left(R_{531}-R_{570}\right)/{\left(R_{731}+R_{570}\right)}^{0.5}$	Canopy	LAI	[43]
RDVI	$(R_{800}-R_{670})/{\left(R_{800}+R_{670}\right)}^{0.5}$	Canopy	Chl, LAI	[44]
REP	$700+[40(R_{670}+R_{780})/2-R_{700}]/\left(R_{740}-R_{700}\right)$	Leaf	Chl	[45]
RVI	$R_{800}/R_{680}$	Canopy	Chl	[31]
SPVI	$1.48\left(R_{800}-R_{670}\right)-1.2\left|R_{530}-R_{670}\right|$	Canopy	Chl, LAI	[25]
TCARI	$3[\left(R_{700}-R_{670}\right)-0.2\left(R_{700}-R_{550}\right)\left(R_{700}/ R_{670}\right)$	Canopy	Chl	[30]
TCARI_ OSAVI	TCARI/OSAVI	Canopy	Chl	[30]
TCARI2	$3[R_{750}-R_{705}-0.2\left(R_{750}-R_{550}\right)\left(R_{750}/ R_{705}\right)$	Canopy	Chl	[46]
TVI	$1.2[1.20\left(R_{750}-R_{550}\right)-2.5\left(R_{670}-R_{550}\right)$	Canopy	LAI, Chl	[47]
VOG1	$R_{740}/R_{720}$	Leaf	Chl, LAI, WC	[48]
VOG2	$\left(R_{734}-R_{747}\right)/\left(R_{715}+R_{726}\right)$	Leaf	Chl, LAI, WC	[48]
VOG3	$\left(R_{734}-R_{747}\right)/\left(R_{715}+R_{720}\right)$	Leaf	Chl, LAI, WC	[48]
WI	$R_{900}/R_{970}$	Soil	WC	[43]
WI/NDVI	$(R_{900}/R_{970})/\left[\left(R_{800}-R_{680}\right)/\left(R_{800}+R_{680}\right)\right]$	Canopy, soil	Chl, LAI, WC	[49]

).

2.5. Modelling

The maximum entropy model (MAXENT), which works on the basis of presence data, was used to compare the suitability of the indices. The presences, in this case, were the centroids of digitized waterlogged areas within all the areas of interest. The environmental conditions were individual hyperspectral indices.

The MAXENT model is often used in a standard parameter setting, which rarely yields the best results [50,51,52], especially in cases where the number of observations is low. Therefore, following the methodology of Warren and Seifert [53] and recommendations [54], firstly, we removed highly correlated variables (r > 0.7) and obtained reduced sets of 9–12 uncorrelated input variables from the original set of 36 variables. Using the Enmevaluate tool, we evaluated all combinations of possible dependencies (linear, quadratic, combined, threshold, stepwise) and selected the one with the lowest AIC (delta AIC = 0). Since the number of waterlogged areas in the areas of interest was small (<50), we used the leave-one-out method for cross-validation. We set the number of replications to the number of waterlogged regions in individual cases. We worked with the original resolution of hyperspectral images 5 m × 5 m when all preparatory work was completed in the ArcGIS 10.8 environment (ESRI 2019. ArcGIS Desktop: Release 10.8. Redlands, CA, USA: Environmental Systems Research Institute).

The entire procedure consists of several steps which, including the processing of results, are shown in Figure 7:

• Finding correlation of input layers and removing highly correlated variables for all territories and each period separately; for our purposes, a value of 0.7 was chosen. Correlation was performed via the freely downloadable ArcGIS SDMToolkit extension.
• A set of detection points (centroids of waterlogged areas) and representative uncorrelated rasters from the previous step were prepared for each territory and each period.
• Using the Enmevaluate tool [53] in the R statistical program environment, variants of the input settings of the MAXENT program were tested precisely according to the procedure described on the website of the Integrative Evolutionary and Conservation Biology Lab (https://sites.google.com/site/thebantalab/). The best input parameters were those with the lowest Akaike information criterion (AIC) value. The result was a correlation parameter (linear, quadratic, combined, threshold, step, or combined) and a regularization parameter.
• The MAXENT application was set according to the results of the previous step and launched.

The Method of Processing Results from the MAXENT Application

The output from the MAXENT program also includes a table evaluating the contributions of individual variables (hyperspectral indices) to the overall probability of occurrence of a waterlogged area. A total of seven cases were processed—four for site 1, each in a different time period, and one for sites 2–4 each, all in the same period.

Let us denote the evaluation of the contributions of individual indices i in a particular case j (combination of time and location) as H_ij. Furthermore, C_ij as a set of indices highly correlated with the index even in the case of j.

For all the k indices that were not selected for modelling in the MAXENT program for the individual case j, we set Hkj as follows:

H_kj = H_ij ∗ cor(k,i),

(1)

where cor represents the correlation of the sought variable k with the selected representative i in the case of j.

The resulting evaluation of the indices was subsequently carried out in several variants:

• All cases, regardless of period and location;
• Spring period, regardless of location;
• Autumn period, regardless of location;
• Site 1 alone, regardless of the period;
• Site 1 alone in spring and autumn;
• The year 2020 alone.

For each variant, an average contribution to the probability of identifying a waterlogged area was obtained.

3. Results

The MAXENT model provides several types of results. The general accuracy of the prediction of waterlogged areas is determined by the value of the area under curve (AUC) or the curve of operating characteristics (ROC). The best in this respect was case 1 in the spring of 2016, followed, in general, by 2017 (

Table 3. Overview of sites, including the setting of binding parameters (L—linear, Q—quadratic, H—threshold) and regularization parameters (rm).

Case	Locality	Period	AUC	Relation	rm
1	1	Spring 2016	0.792 ± 0.235	LQH	2
2	1	Autumn 2016	0.667 ± 0.233	L	3
3	1	Spring 2017	0.715 ± 0.275	LQH	2
4	1	Autumn 2017	0.733 ± 0.307	L	1
5	2	Autumn 2020	0.661 ± 0.258	H	1
6	3	Autumn 2020	0.608 ± 0.291	H	1
7	4	Autumn 2020	0.619 ± 0.337	LQH	1

).

MAXENT also enables the comparison of the significance of variables based on indicators, e.g., a jackknife test of the training set or the test set concerning the gain function (accuracy improvement) or the overall AUC value. We used the expression of the so-called contribution of variability, which shows to what extent individual variables contribute to the general determination of the probability of the observed phenomenon, in this case, the waterlogged area.

Table 4. Percentage of the contributions of hyperspectral indices for identifying waterlogged areas. Bold values represent the values of the variables not processed directly in the given case (period and location) and whose value was derived from a correlated variable according to Equation (1).

Index	Spring 2016	Autumn 2016	Spring 2017	Autumn 2017	Autumn 2020_1	Autumn 2020_2	Autumn 2020_3	Average Contribution
Cari	5.7	55.0	37.8	7.6	14.5	33.7	29.9	26.3 ± 16.6
Dvi	10.1	44.8	30.2	6.6	12.9	29.9	26.2	23.0 ± 12.7
Evi	0.0	0.0	0.0	0.0	37.6	0.0	0.0	5.4 ± 13.2
Mcari	0.4	40.6	26.7	92.4	14.5	25.6	22.3	31.8 ± 27.2
mtvi2	1.1	0.0	0.8	5.7	0.0	62.3	21.1	13.0 ± 21.3
Ndmi	9.2	0.5	1.6	0.0	10.6	0.0	2.2	3.4 ± 4.2
Nmdi	27.7	0.0	6.9	0.0	8.4	2.1	0.0	6.4 ± 9.3
Nvi	54.8	0.0	26.2	6.6	36.4	0.0	0.0	17.7 ± 20.2
Osavi	0.4	0.0	0.0	0.0	0.0	0.0	21.3	3.1 ± 7.4
Pri	0.1	0.1	19.5	5.9	10.6	27.6	25.7	12.8 ± 10.7
Rdvi	0.1	0.0	0.0	0.0	0.0	0.0	21.5	3.1 ± 7.5
tcari2	41.6	3.9	18.8	82.7	0.0	1.9	21.1	24.3 ± 27.4
wi_indvi	0.0	0.0	0.0	0.0	0.0	0.0	46.8	6.7 ± 16.4

shows these contributions in particular periods and locations;

Table 5. Average percentage of the contributions of individual hyperspectral indices for identifying waterlogged areas. Variants: a—all cases, b—spring period, c—autumn, d—site 1, d1—site 1, spring, d2—site 1, autumn, e—the year 2020 (sites 2–4). *** Best, ** second best, * third best.

Index	a	b	c	d	d1	d2	e
cari	** 26.3	* 21.8	** 28.1	* 26.5	* 21.8	*** 66.5	** 26.0
dvi	23.0	20.1	* 24.1	22.9	20.1	0.0	23.0
evi	5.4	0.0	7.5	0.0	0.0	0.0	12.5
mcari	*** 31.8	13.6	*** 39.1	*** 40.0	13.6	* 31.3	20.8
mtvi2	13.0	1.0	17.8	1.9	1.0	0.0	*** 27.8
ndmi	3.4	5.4	2.7	2.8	5.4	3.3	4.3
nmdi	6.4	17.3	2.1	8.7	17.3	3.0	3.5
nvi	17.7	*** 40.5	8.6	21.9	*** 40.5	0.3	12.1
osavi	3.1	0.2	4.3	0.1	0.2	* 43.3	7.1
pri	12.8	9.8	14.0	6.4	9.8	2.9	* 21.3
rdvi	3.1	0.1	4.3	0.0	0.1	0.0	7.2
tcari2	* 24.3	** 30.2	21.9	** 36.7	** 30.2	25.7	7.7
wi_indvi	6.7	0.0	9.4	0.0	0.0	0.0	15.6

then provides the derived aggregated results for the individually processed variants, which are the average of the variable contributions for several cases.

In general, indices relating to chlorophyll content, MCARI (31.8%), CARI (23.6%), and TCARI2 (24.3%), proved to be the best for identifying waterlogged areas. In springtime, the NVI index is the most successful, followed, again, by the indices already mentioned.

Another output of MAXENT is so-called reaction curves, where we can see the effect of a specific value of the monitored variable on the occurrence probability. The reaction curve of the generally best index (according to the analysis)—NVI—shows a significant decrease in probability with values further away from zero. Therefore, waterlogged surfaces are best determined by values with an index close to 0 (see Figure 8). Another example can be the most successful MCARI index from the autumn of 2017, when its contribution to variability was up to 93% (see Figure 9). The threshold here is represented by values greater than 600.

All the MAXENT results can be found in the Figshare Digital Repository (https://doi.org/10.6084/m9.figshare.21710309.v1, accessed on 12 December 2022).

4. Discussion

During the post-war collectivization of the agricultural sector, significant changes occurred in the Czech landscape. The size of blocks of arable land increased, and the average block size is currently greater than 10 ha [55]. This period also saw the homogenization of the landscape caused by the elimination of fragments of woody stands and other nonproductive areas [56,57]. These conditions led to intensified erosion processes and runoff events leading to reservoir silting [58].

The result of landscape homogenization is large blocks of agricultural land, on which waterlogged areas are often present with a minimal or limited production function (Figure 2). At the same time, however, they can be precious from a nonproduction point of view. They retain water in the landscape, contribute to slowing down the runoff, and are centres of increased biodiversity and genetic integrity [8]. They often take on the character of periodic (ephemeral) wetlands. This tends to be a problem when defining their boundaries due to their unstable (periodic/episodic) occurrence [59].

Waterlogging on these sites has various causes. In general, however, it can be said that the leading cause is local physical–geographical factors such as terrain morphology, geological and soil conditions, groundwater level, etc. [60].

The aim of our research was to determine whether it is possible to identify waterlogged areas on agricultural land by means of spectral indices, and with what sensitivity. The analysis made use of spectral indices falling into different categories characterizing the structure, foliage cover, level of chlorophyll and water in plants and in the soil, which could indicate the presence of a waterlogged area. Of course, this range could be expanded to include other spectral indices which are often more difficult to process. However, our intention was to analyse the most commonly used indices, ideally those which can be processed with a certain degree of inaccuracy even without more detailed hyperspectral data using commonly available multispectral data, e.g., from Sentinel-2. The training set was waterlogged areas as they had been identified on the basis of archival and contemporary orthophoto maps.

For our analysis, we used the MAXENT model, which is generally used more for modelling the spatial distribution of species, where it is currently among the most commonly applied models and where it also achieves the most accurate results [52]. This, however, assumes the correct setting of model parameters. For this purpose, we used the Enmevaluate tool, which allows the evaluation of a combination of MAXENT settings based on the lowest Akaike information criterion (AIC). The preselection of input variables for analysis is particularly ambiguous. MAXENT is not susceptible to collinearity; unlike other statistical methods, it is often perceived more as a machine learning tool, and no preselection is necessary [51]. Nevertheless, especially in ecological applications, there is a reduction in input variables, usually through PCA analysis or variable correlation [51], which was also used in our article.

MAXENT allows the creation of models with linear, quadratic, combined, threshold, and step effect of variables on the probability of identification of the observed phenomenon (waterlogged area). The Enmevaluate tool was used as the best setting for input data based on delta AIC depending on the site and period of interest. The type of dependence determines this setting—in our case, it was threshold, quadratic, and linear dependence—and also by choosing the regularization parameter, which is a kind of penalty for overfitting the data. Threshold dependence means a specific threshold value at which linear dependence becomes evident.

Waterlogged areas/water content are/is best detected by three indices considering chlorophyll (MCARI, CARI, and TCARI2). The other indices (DVI, EVI, MTVI2, NDMI, OSAVI, PRI, RDVO, WI_NDVI) were identified as not being very sensitive for the detection of waterlogged areas, i.e., with the lowest average percentage contribution. These less sensitive indices are based on parameters such as plant physiognomy and LAI (for the EVI index).

A high spectral resolution measure, defined as the chlorophyll absorption ratio index (CARI), was developed to minimize the effects of non-photosynthetic materials in the remote estimates of photosynthetically active radiation [61]. Daughtry et al. [29] simplified the CARI index to minimize the combined effects of soil and non-photosynthetic surfaces, leading to the MCARI index.

The modified chlorophyll absorption in reflectance index (MCARI) measures the depth of chlorophyll absorption. It is susceptible to variations in chlorophyll concentration and leaf area index (LAI). MCARI values are not affected by illumination conditions, and the background reflectance from soil and other non-photosynthetic materials was observed. However, Daughtry et al. [29] also showed that the MCARI is influenced by background reflectance and parameters such as LAI and LAI–chlorophyll interaction. Haboudan et al. [30] also showed that the MCARI is sensitive to non-photosynthetic materials, especially at low chlorophyll concentration. That is why the TCARI2 index was created, which compensates for soil reflection and non-photosynthetic components and increases sensitivity at low chlorophyll values.

For the spring period (May), our analysis evaluated the NVI index (new vegetation index) as the most advantageous index. This was initially developed to better distinguish wheat from other crops based on the difference between peak radiances at 777 and 747 nm with reflected radiance at 673 nm in the red area [41]. This index is also described by Liang et al. [62] when studying the leaf area index (LAI) as an indicator of vegetation growth.

The narrow link between relative water content in plants and soil water content had already been revealed by Fischer [63], who suggested plants as an indirect indicator of water content in the soil. Moisture stress inhibits leaf water and chlorophyll content [64,65]. Water stress damages the cell membrane and affects the stability of chlorophyll [66]. Hence, chlorophyll as a parameter may also be used to assess crop water status.

A decreased Chl a/b ratio was also observed under water stress in a study by Zhang et al. [67]. These results indicate that drought can result in damage to the photosynthetic apparatus and may affect other physiological parameters (gas exchange, Chl a/b ratio, water use efficiency, light compensation and light saturation points, maximum photosynthetic rates and dark respiration rates.

The foliar element of the vegetation canopy structure (LAI) aggregates leaf level characteristics and is closely connected with canopy density, physiological processes, and structural traits of leaves for photosynthesis [68]. The value is sensitive to several environmental factors, including moisture [69]. The value increases with moisture, but the relationship to temperature is not straightforward. Global remote sensing also revealed correlations with elevation and precipitation [70,71]. Leaf physiological dynamics are essential indicators of vegetation condition and ecosystem function both from a local and landscape points of view [72,73].

We consider this as the main reason for the indices based on chlorophyll and LAI showing high sensitivity in identifying waterlogged areas. Yet, there are differences between the chlorophyll and LAI-based indices in terms of the detection ability at different times of the year. There is a striking difference between the detection abilities of LAI-based indices in springtime. Chlorophyll-based indices are generally more successful, but achieve the best results in autumn. Both chlorophyll and LAI are important indicators of vegetation conditions and, hence, environmental factors in remote sensing [74,75,76]. However, experimental evaluation of the differences in informative value is beyond the scope of our study, and it is recommended that this be addressed in separate further research.

As part of our analysis, both spectral indices understood as broad-band and narrow-band were used. However, all the indices were derived from hyperspectral data based on the exact wavelengths listed in Table 2. The disadvantage of hyperspectral data is the high purchase price and only partial spatial and temporal coverage of the studied areas.

Some hyperspectral indices can also be derived, to a certain extent, based on multispectral data, e.g., the SENTINEL-2 satellite. Of course, this is only possible if the reflectance in specific wavelengths used for hyperspectral indices can be derived from multispectral data within the considered wavelength bandwidth of the hyperspectral index. In our case, one of the four most advantageous indices—MCARI, CARI, TCARI2, and NVI—can be used with a certain degree of inaccuracy, e.g., the CARI index using wavelengths 700 and 670 nm, which most closely corresponds to band 5—vegetation red edge of the Sentinel-2A satellite (medium wavelength, 704.1 nm, bandwidth, 15 nm), or band 4—red (central wavelength, 664.6 nm, bandwidth, 31 nm). The same is true for the MCARI and tcari2 indices, where, in addition to the wavelengths mentioned for the CARI, a wavelength of 550 nm is also used, which corresponds to band 3—green (central wavelength, 559.8 nm, bandwidth, 36 nm). The NVI index uses wavelengths of 777, 747, and 673 nm, instead of which we can use Sentinel-2A satellite’s band 7—vegetation red edge (central wavelength, 782.8 nm, bandwidth 20 nm), band 6—vegetation red edge (central wavelength, 740.5 nm, bandwidth, 15 nm), and band 4—red (central wavelength, 664.6 nm, bandwidth, 31 nm). With the NVI index, however, this inaccuracy is already quite significant; it is recommended for use only in steps of up to 10 nm [41]; the results, in this case, would be highly debatable.

Despite this, we consider this direction of research to be very promising. In subsequent research, we will focus more closely on the difference between hyperspectral and spectral images.

Of course, the results can also be influenced by precipitation, soil moisture at the time of aerial imaging, and the degree of widespread drought, which can vary considerably in individual years and seasons. Therefore, in addition to spatial diversification, we also processed different periods, covering an above-average rainfall year (2020) and below-average rainfall years (2016 and 2017), as documented (Figure 4). The relative saturation of the soil a few days before the aerial campaign in the spring of 2016 and 2017 was above average, i.e., with sufficient soil moisture, which corresponds to the higher total precipitation of the given month. The autumn periods of 2017 and 2018 were different. While relative soil saturation remained above average in 2016, in 2017, it was the lowest of all the monitored periods and reached around 30%. In terms of precipitation, both years were similar. However, the stormy July of 2016 was reflected in saturation intensity, so the soil did not show any signs of danger from drought. In terms of the accuracy of the MAXENT model forecast, the spring periods were generally more successful (AUC = 0.792 and 0.715, respectively) when soil saturation was generally higher; however, similar success (AUC = 0.733) was also achieved in the autumn of 2017, when the soil suffered from a lack of moisture. Differentiation of the success of individual spectral indices may also relate to this. In the wet spring period, the NVI index sensitive to foliage cover was the most sensitive in terms of waterlogged areas.

On the other hand, indices related to the amount of chlorophyll, distinguishability of bare soil, and soil moisture (OSAVI) performed best in the dry autumn period. The year 2020 was above average in terms of precipitation, which was also evident in the autumn period when the soil did not suffer from a lack of moisture. However, this is precisely why the model results are inconclusive in this period. More extensive aerial campaigns would be needed to confirm the relationship between season, relative soil saturation, and predictive accuracy of models.

5. Conclusions

The overall results can be summarized as follows:

• For the identification of wetland areas, only some of the commonly used hyperspectral indices can be used. In this regard, chlorophyll-based indices worked best for us, with a higher amount of chlorophyll representing a higher probability of a wetland area.
• Although chlorophyll-based indices can be used regardless of the season, better results were achieved in springtime with the NVI index, which represents indices focusing on LAI. A higher probability of occurrence corresponded to low LAI values, indicating low leaf coverage.
• The overall sensitivity of the best indices is statistically significant but does not reach high values. This shows that the use of remote sensing is suitable for the primary selection of wetland areas, which must still be verified in the field.
• The research shows that the exact determination of waterlogged areas in the agricultural landscape is not easy, especially when there is a lack of available hyperspectral data. Such data are not yet readily available in the Czech Republic or other countries. Only recently have the DESIS spectroradiometer on the ISS [77] and two hyperspectral satellites, the Italian PRISMA [78] and the German EnMAP [79], started to supply hyperspectral data. Together with the proposed Copernicus satellite mission CHIME with global coverage [80], we expect hyperspectral data to become more important for various land surface applications. An alternative could be the use of multispectral data using available satellites. However, practical use of such data, due to the lower sensitivity of spectral indices, will be the subject of further research.
• A valuable result, although difficult to interpret biologically, is the role of LAI expressed by the NVI index and the role of chlorophyll in identifying wetland areas. Interpretation cannot be achieved at the level of remote sensing, and further research based on field experiments will be necessary for clarification.

The use of remote sensing can be valuable for both policymakers and practitioners as it can be used to make an initial selection of areas where adaptation and retention measures can be planned and implemented. This remote sensing must be followed by more detailed reconnaissance of the terrain where local knowledge of natural and socio-economic conditions can help in the selection and implementation of the most appropriate measures. These can include the construction of small reservoirs or semiarid polders, as well as the design of nature-based measures, such as the revitalization of rivers, wetlands, or periodic ponds, the aim of which is to improve the state of the landscape affected by climate change and more frequent hydrometeorological extremes.