A Machine Learning Approach for Mapping Chlorophyll Fluorescence at Inland Wetlands

Abstract

Wetlands are a critical component of the landscape for climate mitigation, adaptation, biodiversity, and human health and prosperity. Keeping an eye on wetland vegetation is crucial due to it playing a major role in the planet’s carbon cycle and ecosystem management. By measuring the chlorophyll fluorescence (ChF) emitted by plants, we can get a precise understanding of the current state and photosynthetic activity. In this study, we applied the Extreme Gradient Boost (XGBoost) algorithm to map ChF in the Biebrza Valley, which has a unique ecosystem in Europe for peatlands, as well as highly diversified flora and fauna. Our results revealed the advantages of using a set of classifiers derived from EO Sentinel-2 (S-2) satellite image mosaics to accurately map the spatio-temporal distribution of ChF in a terrestrial landscape. The validation proved that the XGBoost algorithm is quite accurate in estimating ChF with a good determination of 0.71 and least bias of 0.012. The precision of chlorophyll fluorescence measurements is reliant upon determining the optimal S-2 satellite overpass time, which is influenced by the developmental stage of the plants at various points during the growing season. Finally, the model performance results indicated that biophysical factors are characterized by greenness- and leaf-pigment-related spectral indices. However, utilizing vegetation indices based on extended periods of remote sensing data that better capture land phenology features can improve the accuracy of mapping chlorophyll fluorescence.

1. Introduction

Monitoring wetland vegetation is one of the major objectives of the remote sensing of the environment due to its strong relationship with the exchange of atmospheric carbon dioxide (CO2) with land. Numerous environmental scientists have recognized chlorophyll fluorescence as a reliable indicator of a plant’s physiological state, and it is directly connected to photosynthesis [1,2,3]. ChF is the visible red and far-red light emitted by photosynthetic green plant tissues in response to photosynthetically active radiation (PAR) absorption, which occurs between 400 and 700 nm. ChF emission is characterized by two broad peaks that stretch from 650 to 800 nm in the red-edge region, with maxima at 690 nm and 740 nm [4]. Traditionally, the Fv/Fm parameter has been used to measure ChF, which is the most commonly used chlorophyll fluorescence measuring parameter worldwide [5,6]. Fv/Fm is typically measured with a pulse-amplitude modulation fluorometer using various active light sources under light-adapted or dark-adapted conditions [7]. It is a sensitive indicator of a plant’s photosynthetic performance, and lower values may signify stress or photoinhibition, as well as photosynthesis downregulation [8,9,10].

For the past few years, scientists have been able to observe and measure chlorophyll fluorescence in terrestrial ecosystems worldwide using low-resolution remote sensing images. This serves as a proxy for plant photosynthesis [11,12]. Nevertheless, recent studies on conventional ChF retrieval methods in terrestrial ecosystems are challenged by vegetation heterogeneity, seasonal dynamics, and variations in the physical environment among different ecosystems [13,14,15,16]. Algorithms attempting to extract patterns from satellite data are still being developed for retrieving and mapping ChF [17,18,19]. The latest machine learning (ML) approaches in ChF retrieval from the medium-resolution Sentinel and Landsat satellite images have been solely extended for inland and coastal waters [20,21,22,23]. Therefore, it is critical to evaluate novel ML methodologies that enable us to analyze and map the spatio-temporal distributions of plant chlorophyll fluorescence at wetlands frequently using satellite data from the Earth Observation Programme.

Wetland ecosystems are crucial for preserving Europe’s biodiversity and mitigating greenhouse gas levels by acting as carbon repositories. Proper ecosystem management of local water resources requires addressing the issue of wetland restoration. The largest and best preserved of its type in Central and Western Europe is the Biebrza Valley, extending in north-eastern Poland. Unfortunately, this region still faces serious problems in the field of water management and protection of the open fen areas [24,25]. Collecting data on the vegetation and biophysical properties of wetlands using satellite imagery is crucial, as it allows for the management of areas that are difficult or impossible to access for on-site observations. Several studies on the application of remote sensing techniques for environmental monitoring at Biebrza Wetlands have been developed [26,27,28]. Thus, finding the appropriate algorithm to estimate chlorophyll fluorescence accurately is crucial for studying the impact of stress on photosynthesis.

The aim of our study is to determine whether it is possible to use such a machine learning algorithm that is applicable to Sentinel-2 imagery to retrieve ChF data and map this biophysical indicator over wetlands. In situ measurements acquired during the growing season were evaluated and applied for retrieving ChF using a machine learning algorithm utilizing the spectral bands. These latter relations of ChF derived from the field and satellite-imagery-based retrievals were examined for cross-validation. The manuscript follows with the next analyses on (1) how the different spectral bands of Sentinel-2 contribute to the ChF model results, (2) the influence of using a machine learning method, and (3) how this study’s findings can be useful for wetland science and for mapping chlorophyll fluorescence using Earth Observation data.

2. Materials and Methods

2.1. Test Site

The study area is located in north-eastern Poland (Figure 2). Biebrza Wetlands is a unique ecosystem being under RAMSAR Convention for peatlands, as well as its highly diversified flora and fauna, especially reeds and sedges with grasses. The Biebrza National Park was established (in 1993) within the study area as a wetland site of global networking program NATURA 2000. The study area covers a total area of 59,223 ha of which 45% are wetlands. The area is mostly flat, with an elevation of approximately 105 m above sea level. The Biebrza River, which runs through the area, is a naturally occurring lowland river that provides a distinctive reference point for lowland valley mires and river floodplains [29]. There are approximately 70 various types of natural and semi-natural plants documented in vegetation reports [30,31]. To better understand the ecological value of wetlands, we focused on several dominant habitats, such as rush, sedge, sedge moss, and reed communities. According to the European Environment Agency’s (EEA) latest report on wetlands, certain areas of the wetlands still require restoration, particularly those that have been degraded. Our study investigated three main types of plant communities found within the wetland ecosystem: peatlands dominated by sedges, sedge mosses, and reeds (Figure 1). The research site is situated in one of the coolest regions in Poland, where the average daily air temperature throughout the year is only 6.6 °C due to the prevailing temperate continental climate. The coldest month is January, when the average air temperature is around −4.2 °C, and the warmest month is July, with an average temperature of 17.5 °C. Snow coverage usually lasts for a maximum of 140 days, and the average annual precipitation ranges from 550 to 650 mm. Additionally, the growing season in the wetlands is less than 200 days, making it one of the shortest in the country [32].

2.2. Field Measurements

During 2022, five field campaigns were conducted at different times between April and October to collect in situ data, coinciding with Sentinel-2 satellite overpasses around noon local time. Ten sites were selected across various wetland habitats and both sides of the Biebrza River (as shown in Figure 2), and their geographic coordinates were determined using a GNSS receiver. The vegetation’s growth stage was also recorded during each field campaign. The OS5p+ Pulse Modulated Chlorophyll Fluorometer was used to measure ChF in this study [7]. The instrument registers the maximum quantum yield (Fv/Fm) through a dark-adapted test, which is a ratio that represents the maximum potential quantum efficiency of photosystem II (PSII) if all capable reaction centers were open. It estimates the maximum portion of absorbed quanta used in PSII reaction centers [33]. For best results, a dark adaptation time of 20–30 min is recommended [34], and the system provides a series of ten dark-adaption white clips to measure this. After dark adaptation, the end of the fiber optic bundle is placed in the cylindrical opening, and the dark slide of the clip is opened to expose the sample to the fiber optic bundle. Optimal values of Fv/Fm for most land plant species vary within the range of 0.79 to 0.83, while lower values indicate plant stress, particularly photoinhibition [2,8].

2.3. Satellite Data Acquisition

The study area covers two Sentinel-2 granules: 34UFE and 34UEE available from orbits no. 36, 136, and 79 allowing image acquisition every 5–6 days. Satellite images from Sentinel-2A and Sentinel-2B at processing level-2A were retrieved automatically using Google Earth Engine (GEE), i.e., cloud-based platform that provides geospatial data, tools, and computing power for developers to analyze and visualize the world’s satellite imagery and other geospatial data. Users can access and analyze using various programming languages such as Python and JavaScript. It also provides a suite of tools for processing and analyzing the data, including machine learning algorithms for image classification and time-series analysis [35]. For the current study, JavaScript API was used in Earth Engine Code Editor. The average cirrus and cloud cover for the granules in 2022 was 52.15% (Figure 3), whilst the climate in the region is moderate with continental elements, which means that the growing season is short, snow and ice last for a long time, and there is a significant amount of cloud cover. As a result, to create comprehensive images of the area, mosaics were created using images captured by the Sentinel-2 satellite. To combine Sentinel-2 (S-2) data into a mosaic, Google Earth Engine (GEE) uses an algorithm that calculates a high-dimensional weighted geometric median. This approach ensures that the spectral relationships across all the S-2 bands are maintained in the resulting mosaic. Cloudless satellite image mosaics were run spanning around ten days before and after ground measurements. Full list of periods selected for temporal mosaicking of Sentinel-2 images in 2022 over study area is presented in

Table 1. Sentinel-2 imaging periods used for generating cloudless median mosaics.

Name	Start Date	End Date	No. of Sentinel-2 Acquisition Dates Used for Mosaic
April–May	25 April 2022	12 May 2022	7
June–July	20 June 2022	7 July 2022	9
July–August	20 July 2022	4 August 2022	7
August–September	20 August 2022	8 September 2022	10
September–October	20 September 2022	28 October 2022	14

.

2.4. Vegetation Index Calculation

To map the vegetation condition of wetlands using Sentinel-2 data, a total of forty vegetation indices (VIs) were selected and calculated from S2 median mosaics covering the study area. Initially, the most commonly used spectral vegetation indices that describe various aspects of plant growth and health, such as greenness, leaf chlorophyll content, pigments, and water content in the canopy, were identified for the study. These VIs were then arranged into four groups and listed in

Table 2. List of spectral-vegetation-index-dedicated Sentinel-2 bands used for the study.

Application	Abbreviation	Name	Equation	Citation
Spectral Indices of Greenness	AFRI1600	Aerosol Free Vegetation Index 1600	$\mathrm{A}\mathrm{F}\mathrm{R}\mathrm{I}\mathrm{I}1600=R_{842}-0.66\times \frac{R_{1610}}{R_{842}+{0.66\times R}_{1610}}$	[37]
	ARVI	Atmospherically Resistant Vegetation Index	$\mathrm{A}\mathrm{R}\mathrm{V}\mathrm{I}=\frac{R_{842}-R_{665}-0.069\times {(R}_{665}-R_{490})}{R_{842}+R_{665}-0.069\times {(R}_{665}-R_{490})}$	[38]
	CTVI	Corrected Transformed Vegetation Index	$\mathrm{C}\mathrm{T}\mathrm{V}\mathrm{I}=\frac{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}+0.5}{|\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}+0.5|}\times \sqrt{|\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}+0.5|}$	[39]
	EVI	Enhanced Vegetation Index	$\mathrm{E}\mathrm{V}\mathrm{I}=2.5\times \frac{R_{842}-R_{665}}{R_{842}+6.0\times R_{665}-7.5\times R_{490}+1.0}$	[40]
	GDVI	Green Difference Vegetation Index	$\mathrm{GDVI}=R_{842}-R_{560}$	[41]
	GI	Greenness Index	$\mathrm{G}\mathrm{I}=\frac{R_{560}}{R_{665}}$	[42]
	GNDVI	Green Normalized Difference Vegetation Index	$\mathrm{G}\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\frac{R_{842}-R_{560}}{R_{842}+R_{560}}$	[43]
	mNDVI	Modified NDVI	$\mathrm{m}\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\frac{R_{842}-R_{665}}{R_{842}+R_{665}-2.0\times R_{443}}$	[44]
	NDVI	Normalized Difference Vegetation Index	$\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\frac{R_{842}-R_{665}}{R_{842}+R_{665}}$	[45]
	rNDVI	Renormalized Difference Vegetation Index	$\mathrm{r}\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\frac{R_{842}-R_{665}}{\sqrt{R_{842}+R_{665}}}$	[46]
	NDRE	Normalized Difference NIR / Red Edge	$\mathrm{N}\mathrm{D}\mathrm{R}\mathrm{E}=\frac{R_{842}-R_{705}}{R_{842}+R_{705}}$	[47]
	PPR	Normalized Difference 550/450 Plant Pigment Ratio	$\mathrm{P}\mathrm{P}\mathrm{R}=\frac{R_{560}-R_{443}}{R_{560}+R_{443}}$	[48]
	PVR	Normalized Difference 550/650 Photosynthetic Vigour Ratio	$\mathrm{P}\mathrm{V}\mathrm{R}=\frac{R_{560}-R_{665}}{R_{560}+R_{665}}$	[48]
	RENDVI	Red Edge Normalized Difference Vegetation Index	$\mathrm{R}\mathrm{E}\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\frac{R_{740}-R_{705}}{R_{740}+R_{705}}$	[49,50]
	SAVI	Soil Adjusted Vegetation Index	$\mathrm{S}\mathrm{A}\mathrm{V}\mathrm{I}=\frac{R_{842}-R_{665}}{R_{842}+R_{665}+0.725}\times (1.0+0.725)$	[51]
	SLAVI	Specific Leaf Area Vegetation Index	$\mathrm{S}\mathrm{L}\mathrm{A}\mathrm{V}\mathrm{I}=\frac{R_{842}}{R_{665}+R_{2190}}$	[52]
	SR	Simple Ratio 842/665	$\mathrm{S}\mathrm{R}=\frac{R_{842}}{R_{665}}$	[53,54]
	SRT	Simple Ratio 1610/2190	$\mathrm{S}\mathrm{R}\mathrm{T}=\frac{R_{1610}}{R_{2190}}$	[55]
	S2REP	Sentinel-2 Red-Edge Position Index	$\mathrm{S}2\mathrm{R}\mathrm{E}\mathrm{P}=705+35\times \frac{(\left(R_{665}+R_{783}\right)÷2)-R_{705}}{(R_{740}-R_{705})}$	[56]
Leaf Chlorophyll Content	CCCI	Canopy Chlorophyll Content Index	$\mathrm{C}\mathrm{C}\mathrm{C}\mathrm{I}=\frac{{(R}_{842}-R_{705})÷{(R}_{842}+R_{705})}{{(R}_{842}-R_{665})÷{(R}_{842}+R_{665})}$	[57]
	CVI	Red-edge-band Chlorophyll Index	$\mathrm{C}\mathrm{V}\mathrm{I}=\frac{R_{842}\times R_{664}}{\sqrt{R_{560}}}$	[58]
	IRECI	Inverted Red-edge Chlorophyll Index	$\mathrm{I}\mathrm{R}\mathrm{E}\mathrm{C}\mathrm{I}=\frac{R_{842}-R_{665}}{R_{740}÷R_{705}}$	[56]
	LCI	Leaf Chlorophyll Index	$\mathrm{L}\mathrm{C}\mathrm{I}=\frac{R_{842}-R_{705}}{R_{842}+R_{665}}$	[59]
	MCARI	Modified Chlorophyll Absorption in Reflectance Index	$\mathrm{M}\mathrm{C}\mathrm{A}\mathrm{R}\mathrm{I}=1.2\times \left[2.5\times \left(R_{842}-R_{665}\right)-1.3\times \left(R_{842}-R_{560}\right)\right]$	[60]
	TCARI	Transformed Chlorophyll Absorption Ratio	$\mathrm{T}\mathrm{C}\mathrm{A}\mathrm{R}\mathrm{I}=3\times \left[\left(R_{705}-R_{665}\right)-0.2\times \left(R_{705}-R_{560}\right)\right]\times \left(R_{705}÷R_{665}\right)$	[61]
	TCI	Triangular Chlorophyll Index	$\mathrm{T}\mathrm{C}\mathrm{I}={1.2\times (R}_{705}-R_{560})-1.5\times {(R}_{665}-R_{560})\times \sqrt{{(R}_{705}÷R_{665})}$	[62]
Leaf Pigments	ARI	Anthocyanin Reflectance Index	$\mathrm{A}\mathrm{R}\mathrm{I}=\frac{1}{R_{560}}-\frac{1}{R_{705}}$	[63]
	BGI	Blue Green Pigment Index	$\mathrm{B}\mathrm{G}\mathrm{I}=\frac{R_{443}}{R_{560}}$	[64]
	BRI	Browning Reflectance Index	$\mathrm{B}\mathrm{R}\mathrm{I}=\frac{\frac{1}{R_{560}}-\frac{1}{R_{705}}}{R_{842}}$	[65]
	CI	Coloration Index	$\mathrm{C}\mathrm{I}=\frac{R_{665}-R_{490}}{R_{665}}$	[58]
	GLI	Green Leaf Index	$\mathrm{G}\mathrm{L}\mathrm{I}=\frac{{2\times R}_{560}-R_{665}-R_{490}}{{2\times R}_{560}+R_{665}+R_{490}}$	[66]
	PBI	Plant Biochemical Index	$\mathrm{P}\mathrm{B}\mathrm{I}=\frac{R_{842}}{R_{560}}$	[67]
	PSRI	Plant Senescence Reflectance Index	$\mathrm{P}\mathrm{S}\mathrm{R}\mathrm{I}=\frac{R_{665}-R_{490}}{R_{720}}$	[68]
	SIPI	Structure Insensitive Pigment Index	$\mathrm{S}\mathrm{I}\mathrm{P}\mathrm{I}=\frac{R_{842}-R_{443}}{R_{842}+R_{665}}$	[69]
Canopy Water Content	GVMI	Global Vegetation Moisture Index	$\mathrm{G}\mathrm{V}\mathrm{M}\mathrm{I}=\frac{{(R}_{842}+0.1)-{(R}_{2190}+0.2)}{{(R}_{842}+0.1)+{(R}_{2190}+0.2)}$	[70]
	MSI	Moisture Stress Index	$\mathrm{M}\mathrm{S}\mathrm{I}=\frac{R_{1610}}{R_{842}}$	[71]
	NDII	Normalized Difference Infrared Index	$\mathrm{N}\mathrm{D}\mathrm{I}\mathrm{I}=\frac{R_{842}-R_{1610}}{R_{842}+R_{1610}}$	[72]
	NDWI	Normalized Difference Water Index	$\mathrm{N}\mathrm{D}\mathrm{W}\mathrm{I}2=\frac{R_{842}-R_{2190}}{R_{842}+R_{2190}}$	[73]
	NMDI	Normalized Multi-band Drought Index	$\mathrm{N}\mathrm{M}\mathrm{D}\mathrm{I}=\frac{R_{865}-(R_{1610}-R_{2190})}{R_{865}+{(R}_{1610}-R_{2190})}$	[74]
	SIWSI	Shortwave Infrared Water Stress Index	$\mathrm{S}\mathrm{I}\mathrm{W}\mathrm{S}\mathrm{I}=\frac{R_{865}-R_{1610}}{R_{865}+R_{1610}}$	[75]

.

2.5. Model Establishment and Evaluation

The XGBoost algorithm was investigated to build a model for estimation of chlorophyll fluorescence, which is an open-source software library for gradient boosting algorithms that is designed to be highly efficient, flexible, and scalable. It is a popular machine learning framework that is widely used in industry and academia for supervised learning tasks, such as classification and regression. The XGBoost algorithm works by combining multiple weak models to create a strong model. It trains each model iteratively by adjusting the weights of misclassified samples, with the aim of minimizing a loss function. The final model is a weighted sum of all the models, with the weights determined by their performance on the training data [76]. XGBoost is known for its speed and performance, and it supports regularization, uses the power of parallel processing, handles missing values, and performs cross-validation [77].

The input for the model was an array containing r rows and c columns, where r represents the total number of observations with available data on predictors (vegetation indices) and references (Fv/Fm), and c represents the number of columns with predictors calculated for each growth stage. All the predictors were scaled to a range of values between zero and one. The proposed model in supervised learning referred to the mathematical structure by which the prediction was made from the satellite-based vegetation indices (listed in Table 2). A linear prediction model was created, which calculates predictions based on a combination of weighted satellite-derived data and parameters. The available dataset was split into a training set and a validation set in a ratio of 80:20. To improve the accuracy of the prediction model, the XGBoost algorithm was fine-tuned using tenfold cross-validation and learning curve methods during the training phase to determine the best parameter configuration.

In order to evaluate the performance of the model, a cross-validation was carried out using the leave-one-out approach, which is a variant of k-fold cross-validation, where k is equal to the number of observations in the dataset. It should be noted that the feature selection process was repeated at each iteration to prevent the model from benefiting from knowledge of the data used for validation. The model performance was assessed using various metrics, including the coefficient of determination (R²), which is the square of the correlation coefficient, as well as the mean absolute error (MBE), root mean square error (RMSE), and relative root mean square error (RRMSE).

3. Results

3.1. Analysis of Ground-Measured Chlorophyll Fluorescence

The variations in maximum quantum yield (Fv/Fm) characterizing chlorophyll fluorescence in different vegetation growth periods at the wetlands in 2022 are shown in Figure 4. The most ChF records were made in late June and early July (115), while due to unfavorable weather conditions in the period September–October, only 48 field measurements were carried out. The Fv/Fm reached its maximum value (0.82) at the late flowering stage in August–September. It was noted at its lowest at the late beginning and at the end of the growing season (0.76), in the periods April–May and September–October, respectively. The lowest standard deviations appeared in the mid-season (July–August–September), at 0.010–0.011, whilst the highest standard deviations appeared in the early (April–May) and late stages (September–October) of wetland vegetation development, at 0.026 and 0.029, respectively.

3.2. Features Selected by Algorithm for Mapping ChF

The initial feature set is composed of 19 greenness-related vegetation indices, 7 indices characterizing leaf chlorophyll content, 8 features concerning leaf pigments, and 6 indices on canopy water content. Each feature parameter in the feature set is analyzed individually, and Figure 5 shows the correlation coefficient between each of these parameters. It can be noted that the overall features have high positive and negative correlations with the other ones. The water content, leaf pigments, and canopy structure of herbaceous plants in wetlands exhibit clear changes that can be used to monitor their growth dynamics. Additionally, measuring the chlorophyll content of the canopy is a sensitive method for accurately determining vegetation growth patterns. Among all analyzed features, the indices based on Sentinel-2 spectral bands of near-infrared (NIR) and shortwave-infrared (SWIR) received the most correlations among each other. The parameters related to the water content of the canopy, such as the area and amplitude of the edges in the near-infrared (NIR) and shortwave-infrared (SWIR) regions, exhibit significant changes in response to variations in vegetation conditions caused by river water-level fluctuations and spring floods. The algorithm selects the combination of wavelenghts with the least collinearity, and then the correlation coefficient between each spectral index based on combinations of various wavelengths with another one shows non-uniformity.

In this study, we examined the use of XGBoost for automating the feature selection process on the initial feature set. The outcome of this process is presented in Figure 6. It can be observed that the Simple Ratio on the 1610 and 2190 bands (SRT) has the highest importance, while the importance of the ARVI is slightly lower. Both indices reflect vegetation growth dynamics and temporal changes in photosynthetic capacity. The vegetation indices ARI, CI, and BRI, which reflect plant senescence, can detect even minor changes in leaf pigments and capture the phenomenon of vegetation browning. When herbaceous plants are exposed to various environmental stresses such as seasonal flooding, high variations in salt levels, low soil oxygen levels, and waves over time, their cell viability and biochemical composition change, which can affect the morphology of leaves, the distribution of leaf inclination, and the structure of the canopy. Therefore, the SRT, which shows healthy and unhealthy vegetation, has a distinct response to environmental conditions. According to a previous study [34], it was found that under saturation conditions, the variation in far-red fluorescence is mainly influenced by the structure of the canopy. The variability in the broadband light spectrum ranging from 641 to 800 nm, which is included in the ARVI and ARI indices, is mainly influenced by the optical properties of leaves and the structural characteristics of the canopy [15]. Finally, these spectral vegetation indices can reflect the vigor of wetland vegetation linked to the changes in leaf and canopy structure and subsequently to chlorophyll fluorescence.

3.3. Evaluation of the Model Performance

To evaluate the XGBoost algorithm results, we compared the predicted chlorophyll fluorescence Fv/Fm with reference data, collected in situ during the growing season in 2022. The evaluation of the model performance is presented in Figure 7. It is worth noting that high coefficients of determination at R² = 0.71 and low errors at the levels of RMSE = 0.012, rRMSE = 0.016, and MAE = 0.009 were found with the used algorithm. These findings indicate quite good prediction, especially having the short time span of analysis (April–October 2022). The evaluation results indicate that the XGBoost algorithm has high monitoring accuracy, which is consistent with the results of a previous study [78] that used different models. This suggests that the XGBoost algorithm is not only universal but also scalable for vegetation monitoring.

3.4. Spatial and Temporal Patterns of Chlorophyll Fluorescence

The estimated spatial changes in ChF from April to October 2022 are shown in Figure 8. The areas with non-wetland vegetation such as forests, arables, and meadows were excluded using classification from previous studies [79,80]. Most of the greening areas are in the Biebrza Lower Basin, consistent with the previous studies [26,27,32]. In addition, the spatial estimations on wetland vegetation with high ChF also confirm the results on mapping gross primary productivity (GPP) with remote sensing data from 2022 [80].

Spatio-temporal changes for five periods are presented in Figure 9. The periods of April–May and June–July exhibit a significant increase in vegetation cover in the Biebrza Valley, particularly in the wild regions of the Biebrza Lower Basin, whilst at the turn of the months July–August, the majority of the Biebrza Valley shows insignificant ChF decreases, with some noticeable browning. This is in agreement with the study of Okruszko, 1990 [29]. During the following month August until October, areas with significant ChF changes were dominated by browning trends under the temperate climate. This is also supported by other studies on Cepkeliai wetlands in Lithuania that are located 200 km from the Biebrza Valley [81]. Moreover, Simanauskiene et al. 2019 [81] found that the summit of the growing season of herbaceous plants is noted in June, while vegetation index variations in July reflect environmental changes indicating browning. The findings suggest that there are significant fluctuations in the spatial distribution of ChF changes across different time periods. We observe an overall greening development stage at a maximum ChF of 0.84 and browning at 0.70.

4. Discussion

The research on a time series of combined optical Sentinel-2 and field remote sensing data led to providing maps of chlorophyll fluorescence over wetlands in a large area. So far, several studies have documented the potential of remote sensing data to map chlorophyll fluorescence [13,16,79]. Whilst Sinha et al. [18] used medium-resolution remote sensing imagery for estimating seasonal fluorescence dynamics at deciduous forests, we mapped the areas difficult to access with high water-level fluctuations. Therefore most studies examined relatively homogeneous areas such as forests; however, larger heterogeneous landscapes with high biodiversity and conservation practices pose a challenge. Nevertheless, the application of Sentinel-2 imagery to estimate chlorophyll fluorescence at wetlands remains underexplored. In this context, we endeavor to figure out the potential of Sentinel-2 multispectral data for estimating the plant state reflecting seasonal variations in phenology at the areas being protected under the Ramsar Convention in Poland. We studied the XGBoost algorithm supported by the combination of Fv/Fm field data to select the most appropriate satellite-derived features for the spatio-temporal monitoring of chlorophyll fluorescence.

Our study demonstrates the potential of 10 m satellite imagery for estimating the chlorophyll fluorescence parameter, which could not be assessed directly from remote sensing data. Comparing coarse-resolution MODIS satellite imagery, Simanauskiene et al. 2019 [81] found that commonly used NDVI values can be recognized as an appropriate indicator for finding degraded peatland areas. However, early detection of degrading peatlands with high spatial and temporal resolution still meets the challenge. Moreover, a suitable algorithm for accurate chlorophyll fluorescence estimation to find the effect of stresses on the photosynthetic process is still required for management practices [82,83,84]. Our study showed a novel approach for mapping ChF and subsequently indicated spatio-temporal greening and browning variations. The features selected by XGBoost enable us to identify plant phenology shifts, such as greening and browning.

Many studies have used the XGBoost algorithm for monitoring environmental changes with remote sensing data [22,23,85]. Jing et al. 2022 [85] found that data collected from one year provide satisfactory results; however, in order to improve the accuracy of prediction on features by XGBoost, more data characterizing meteorological and environmental variables could be included. In our study, spectral vegetation indices and field data were also collected as features with a higher importance in one analyzed year. We are aware that the more field data we collect in subsequent years, the higher accuracy of the model we might assess. The reasons we applied XGBoost are as follows: small dataset (one-year collection), model architecture (decision trees), structured data (Fv/Fm and vegetation indices), and feature importance scores (straightforward and easy to understand). XGBoost and deep learning models are both machine learning algorithms, but they differ in several ways: (a) Model architecture: XGBoost uses an ensemble of decision trees as base learners, whereas deep learning models use artificial neural networks that are composed of multiple layers. (b) Input data format: XGBoost is well-suited for structured data that are arranged in rows and columns, while deep learning models can handle both structured and unstructured data, such as images and text. (c) Computational requirements: deep learning models typically require more computational resources than XGBoost, such as GPUs or TPUs, and can take longer to train. (d) Interpretability: XGBoost provides more interpretability than deep learning models, as it can output feature importance scores and decision rules. In contrast, deep learning models are often considered “black box” models, as it can be difficult to interpret them and understand how they arrive at a particular prediction. (e) Performance on small datasets: XGBoost can perform well on small datasets, whereas deep learning models typically require large amounts of data to perform well.

We utilized Sentinel-2 median mosaics due to cloud coverage frequency at the study site. This process allowed us to reduce the dimensionality of the input array consisting of satellite daily images and adjust the periods of temporal mosaics covering field campaigns at the Biebrza Valley. Results on the quality of the compositions from Sentinel-2 satellite images by Shepherd et al. 2020 [86] confirmed the use of improved cloud-free and composed daily, weekly, or monthly mosaics for regular land monitoring. Regarding unfavorable weather conditions in the study area, Sentinel-2 median mosaics were taken into account as well. In our study, we did not explicitly consider the impact of view zenith angle difference and its derivative consequences on the model performance in temporal pattern prediction [87]. However, this is an important factor to consider in future studies. Therefore, we recommend that future studies investigating temporal patterns in remote sensing data should consider this factor and explicitly evaluate its impact on model performance.

However, some limitations must be mentioned in next future studies. Even though the XGBoost algorithm effectively detected ChF patterns at the local site and the method presented might be applied to other areas, the parameters require further optimization. In this context, the XGBoost model will be developed applying additional parameters characterizing high water-level fluctuations, e.g., surface roughness and soil moisture establishments from satellite radar data. Additionally, regarding weather and optical satellite daily image constraints, the obtained median mosaics were few, and the coverage was limited. Therefore, this study did not focus on large-scale environmental monitoring or detecting stress periods. However, including meteorological information such as air temperature, precipitation, and humidity derived from ERA-5 Land reanalysis could be valuable for future research.

Our study aimed to emphasize the possibility of applying free-of-charge Sentinel-2 imagery for mapping chlorophyll fluorescence with fine resolution. Considering the reports on the state of European wetlands from the EEA and the anthropogenic influence on the wetland biodiversity and sustainable management of the Biebrza Valley in Poland [25,84], the proposed study might support biodiversity conservation and management practices that are the subject of investigations from other academic units [88,89].

5. Conclusions

Our study demonstrates the feasibility of mapping chlorophyll fluorescence for an area of terrestrial wetlands using remote sensing indices derived from the Sentinel-2 satellite with high accuracy (R² = 0.71, RMSE = 0.012). However, fluorescence information has so far most often been related point-wise to specific study locations. We have shown that by using ground-based field measurements, it is possible to use machine learning techniques (XGBoost algorithm) to relate the solution spatially, enabling the analysis of wetland conditions and wetland phenological changes (growth, full vegetation development, dieback). The conducted feature selection showed that among the most important remote sensing indices in the wetlands were a group of leaf pigment (Anthocyanin Reflectance Index, Coloration Index, Browning Reflectance Index, Blue Green Pigment Index) and water indices (Normalized Multi-band Drought Index, Shortwave Infrared Water Stress Index, Moisture Stress Index). In detail, the group of water indices shows a strong relationship between water content and chlorophyll fluorescence, which demonstrates the high dependence of photosynthesis on water conditions in wetlands, which are an important absorber of atmospheric carbon dioxide, and drainage could adversely affect this potential. Generally, regarding our results and F1 scores presented in Figure 6, the Simple Ratio Index based on Sentinel-2 Short Wave Infrared (B11 and B12) bands revealed the highest importance in the model. AFRI1600 based on Infrared and Short Wave Infrared (Sentinel-2 B8 and B11) is among the crucial features to model as well. However, we are aware that other vegetation indices, such as the ARVI and ARI, which are based on spectral reflections in the visible blue and green bands, also significantly influenced the estimation of chlorophyll fluorescence. Therefore, it is recommended to conduct analyses using a wide range of indices and other biophysical parameters for research purposes. It is important to note that the specific bands used for fluorescence modeling may vary depending on the type of plant or vegetation being studied, as well as the environmental conditions being assessed. This research is a contribution to further fluorescence mapping studies, especially in the context of the upcoming European Space Agency mission, the FLuorescence EXplorer (FLEX), which will aim to determine fluorescence data for the globe on a continuous basis. Further surveys will allow the proposed solution to be extended to other plant ecosystems (forests, agricultural fields), which is particularly important in view of climate change and sustainable agriculture. A possible future direction for the research will be to combine the fluorescence data with other biophysical variables as well as climatic conditions, identifying the relationship between them and their influence on chlorophyll fluorescence.