Random Forest-Based Soil Moisture Estimation Using Sentinel-2, Landsat-8/9, and UAV-Based Hyperspectral Data

Abstract

Accurate spatiotemporal monitoring and modeling of soil moisture (SM) is of paramount importance for various applications ranging from food production to climate change adaptation. This study deals with modeling SM with the random forest (RF) algorithm using datasets comprising multispectral data from Sentinel-2, Landsat-8/9, and hyperspectral data from the CoSpectroCam sensor (CSC, licensed to AgriWatch BV, Enschede, The Netherlands) mounted on an unmanned aerial vehicle (UAV) in Iran. The model included nine bands from Landsat-8/9, 11 bands from Sentinel-2, and 1252 bands from the CSC (covering the wavelength range between 420 and 850 nm). The relative feature importance and band sensitivity to SM variations were analyzed. In addition, four indices, including the perpendicular index (PI), ratio index (RI), difference index (DI), and normalized difference index (NDI) were calculated from the different bands of the datasets, and their sensitivity to SM was evaluated. The results showed that the PI exhibited the highest sensitivity to SM changes in all datasets among the four indices considered. Comparisons of the performance of the datasets in SM estimation emphasized the superior performance of the UAV hyperspectral data (R² = 0.87), while the Sentinel-2 and Landsat-8/9 data showed lower accuracy (R² = 0.49 and 0.66, respectively). The robust performance of the CSC data is likely due to its superior spatial and spectral resolution as well as the application of preprocessing techniques such as noise reduction and smoothing filters. The lower accuracy of the multispectral data from Sentinel-2 and Landsat-8/9 can also be attributed to their relatively coarse spatial resolution compared to the CSC, which leads to pixel non-uniformities and impurities. Therefore, employing the CSC on a UAV proves to be a valuable technology, providing an effective link between satellite observations and ground measurements.

1. Introduction

Soil moisture monitoring is crucial in various fields, including soil science [1], ecology [2], hydrology [3,4], and water resource management [5]. In agriculture, accurate monitoring of SM at different stages of crop growth is particularly important for water stress detection [6] and crop yield estimation [7]. In addition, SM is considered an important component of the hydrological cycle as it effectively controls the ratio of rainfall conversion to infiltration and runoff [8,9]. Consequently, spatial and temporal monitoring of SM changes can be an important prerequisite for climate change mitigation and adaptation [10], which is essential not only for the sustainability of agricultural systems, but also for the further development of precision agriculture and for ensuring food security [11,12,13,14,15].

Direct methods for estimating SM are limited to point-based approaches that provide localized SM measurements at specific depths and points. However, such point-based measurements do not represent distribution patterns of SM, as SM is highly variable both spatially and temporally [16,17]. Point-based methods also have disadvantages, such as the high cost of equipment, the high labor input for field-level measurements, and the time-intensive nature of soil sampling. Therefore, an alternative solution for SM estimation is the use of remote sensing technology [18,19,20]. This approach offers a cost-effective way to achieve both temporal and spatial coverage, e.g., for precision agriculture, and to avoid the above-mentioned problems of point-based soil sampling.

In general, an increase in SM leads to a decrease in soil reflectance [21]. This complex relationship can lead to significant errors in modeling and remote sensing-based estimates of SM. Numerous studies have addressed the modeling of SM using multispectral satellite data. A study from South Africa highlighted the high potential of Sentinel-1 and Sentinel-2 data for SM estimation [22]. The findings indicated that RF outperformed support vector regression (SVR) and simple linear regression algorithms in estimating SM. In another study, SM estimation from Sentinel-2 data was performed using a convolutional neural network (CNN) in Austria and Australia [23]. The results showed that the short-wave infrared-1 band (SWIR-1) had the highest sensitivity, while the blue band had the lowest sensitivity for SM. The significant potential of Landsat-8 data for SM estimation has been demonstrated in various studies [24,25,26]. Despite the considerable potential of multispectral satellite imagery for SM estimation, its limited spatial resolution hampers its applicability at smaller scales. The increasing use of UAVs has become prominent due to their enhanced temporal flexibility and higher spatial resolution compared with satellite data. Recent advances include the use of small digital sensors on UAVs that provide exceptional spectral and spatial resolution [27]. There is a growing body of research on SM prediction through the integration of UAVs and various sensors. Notable studies have used multispectral sensors on UAVs [28,29,30], while a few have integrated hyperspectral sensors [31,32]. Modeling SM from these datasets requires the extraction of relationships between SM content and surface reflectance and emission. To achieve this, machine learning algorithms prove helpful in establishing the intricate connections within the data.

The versatility of machine learning algorithms has been demonstrated on a variety of remote sensing data, including optical, radar, and UAV data [33]. Comparative analyses with other parametric approaches, including iterative and Bayesian methods, consistently show the superior performance of machine learning algorithms in numerous studies [34]. In contrast to linear mathematical models, machine learning is characterized by the detection of non-linear relationships between predictor and target variables, as has been observed in various studies [31]. A distinct advantage of machine learning methods is that they rely less on preprocessing data for modeling, especially due to their ability to process different datasets without extensive transformations [35,36]. The machine learning landscape includes several established methods such as the artificial neural network (ANN), RF, SVR, extreme learning machine (ELM), partial least squares regression (PLSR), and Gaussian process regression (GPR). Numerous studies have applied machine learning algorithms to estimate SM using various remote sensing data and demonstrated their effectiveness [24,31,37,38,39]. In a study from China, RF, KNN, PLSR, and backpropagation neural network (BPNN) models were used to estimate SM in a maize field using multispectral, RGB, and thermal data collected by a UAV-mounted sensor [40]. The results showed that multispectral data consistently provided more accurate results than RGB and thermal data, regardless of the model used. Furthermore, the results emphasized the superior performance of the RF model in SM estimation compared to KNN, PLSR, and BPNN. In another study, UAV-based RGB data were used to estimate SM in a bare soil field in Iran using RF, M5 tree (M5P), multilayer perceptron (MLP), and sequential minimal optimization regression (SMOreg) algorithms [41]. The results showed that the RF model was superior to the other models in accurately estimating SM. In addition, further research addressed SM modeling using the general linear model (GLM) at depths of 2, 8, and 12 inches in South Georgia using UAV hyperspectral data [42]. The results showed a statistically significant correlation between hyperspectral data and ground data at a depth of 2 inches. However, the models tested for depths of 8 and 12 inches were not able to estimate SM with the same accuracy. Similarly, another study focused on developing a quantitative model to estimate SM using ELM and RF models [31]. A hyperspectral sensor mounted on the Matrice 600 UAV captured electromagnetic waves in the wavelength range from 400 to 1000 nm. The comparison of the two machine learning algorithms showed the superiority of the RF model over the ELM. Estimating SM from UAV thermal and multispectral sensors, Bertalan [43] used the RF, GLM, elastic net regression (ENR), and robust linear model (RLM) models. The results showed that the RF model, utilizing both thermal and multispectral data, outperformed the other models. Comparing the accuracy of machine learning models in estimating SM, de Oliveira [44] showed that the RF model outperformed the average neural network and weighted k-nearest neighbor (KNN) models in spatial and temporal modeling of SM. While several studies have focused on SM modeling using satellite or UAV data, few comparative analyses directly evaluate the accuracy and performance of UAV-based hyperspectral sensors versus satellite-based multispectral data in SM modeling. This research gap prevents a comprehensive understanding of the strengths and limitations of each data source in accurately estimating SM under different environmental conditions and at different scales. Such comparative analyses are critical for selecting the most appropriate data source for specific SM modeling applications, considering factors such as spatial resolution, spectral range, and temporal frequency.

The main objectives of this study are (1) to compare the performance of multispectral data from Sentinel-2 and Landsat-8/9 with hyperspectral data from the UAV-mounted CSC in SM estimation using the RF algorithm, (2) to investigate the relationships between all possible two-band combinations of each dataset with SM, and (3) to determine the most sensitive bands of each dataset in SM estimation.

2. Materials and Methods

2.1. The Study Area

The study was conducted on a research farm of the University of Tehran, located in Pakdasht of Iran (latitude 35°28′25″N, longitude 51°33′35″E), with an area of 100 hectares. This area is slightly sloped and has an average elevation of 985 m above sea level. The region receives an average of 217 mm annual precipitation. The predominant land uses in the study area include bare soil, alfalfa, wheat, and corn. This land use composition can vary from year to year, depending on the availability of water resources. Figure 1 shows an overview of the study area and the location of the soil samples collected for the individual land uses.

2.2. Required Data

2.2.1. In Situ Soil Moisture Data

In each field campaign, 45 points (12 samples in alfalfa, 11 samples in wheat, 11 samples in corn, 11 samples in bare soil) were selected for in situ soil sampling. At each sampling point, five topsoil cores (0–10 cm depth) were combined to create composite samples (45 × 5 soil samples were taken from 45 points in each field campaign). Soil moisture content was measured using the gravimetric method [4].

2.2.2. Soil Properties Data

The spectral reflectance of the soil surface is influenced by various soil properties [45,46,47]. Therefore, 25 points were selected to investigate the spatial variations of electrical conductivity (EC), pH, soil texture, and soil color in the region. The hydrometric method proposed by Bouyoucos [48] was used to assess soil texture at the selected points. In addition, an EC meter was used to measure the EC value, a pH meter to determine the pH of the soil, and the Munsell color system to characterize the soil color. Notably, all samples showed remarkable homogeneity in terms of texture, with silty loam identified as the predominant soil texture in the field. In addition, the samples were uniform in color, pH, and salinity. The premise of this study is therefore that soil properties in the region are predominantly homogeneous, with the spectral response of the soil being primarily influenced by variations in SM.

2.2.3. Satellite Data

The Sentinel-2 Level-2A (L2A) surface reflectance images were obtained from the European Space Agency (ESA) website (https://scihub.copernicus.eu, accessed on 28 March 2017). The Level-2A product provides atmospherically corrected surface reflectance images derived from the corresponding Level-1C products. The Aerosol, Blue, Green, Red, Red Edge 1, Red Edge 2, Red Edge 3, NIR, Red Edge 4, SWIR-1, and SWIR-2 bands were selected for SM modeling as they are particularly sensitive to SM [23].

The specifications of the Landsat-8 and Landsat-9 satellites are similar; the only relevant difference is that the radiometric resolution of 12 bits for Landsat-8 has been upgraded to 14 bits for Landsat-9, and the overall signal-to-noise ratio has been slightly improved. Corrected Collection-2 Level-2 images of these two satellites were downloaded from the United States Geological Survey (USGS) website. The Aerosol, Blue, Green, Red, NIR, SWIR-1, SWIR-2, TIR-1, and TIR-2 bands were used for modeling SM.

2.2.4. UAV Hyperspectral Data Acquisition

In this study, data from a USB2000PLUS hyperspectral sensor (CSC) was used in addition to the satellite data. The CSC has spectroscopic capabilities in the wavelength range from 339.6 to 1028.8 nm, covering the visible and near-infrared (NIR) regions, and offers a spectral resolution of 0.35 nm (with a total of 2048 bands). The DJI Matrice 100 UAV served as the carrier for the CSC. Data acquisition with the CSC was performed on seven different days, synchronized with satellite passes over the study area. To mitigate atmospheric effects, data acquisition was done on sunny days. Figure 2 shows the Matrice 100 UAV and the CSC mounted on it. The availability of datasets in each data collection is shown in Figure 3.

CoSpectrocam software served as an intermediary to establish a connection between the UAV and the CSC. Before the flight, the CSC was calibrated using a Spectralon White Diffuse Reflectance Standard. Equation (1) was used to measure the reflectance of the objects [49]:

$$\rho =\frac{{\mathrm{L}}_{sam}}{{\mathrm{L}}_{ref}}\times {\rho }_{ref}$$

(1)

where ρ is the reflectance of the sample, L_sam and L_ref are the measured radiance of the sample and the white reference, respectively, and ρ_ref is the reflectance of the white standard reference.

DJI Go software (version 3.1.61) was used to establish a connection between the UAV and the remote control. Drone Harmony software (version 2.3.0, Drone Harmony Company, Zurich, Switzerland) was used for the automatic navigation of the UAV. The flight parameters were set using Drone Harmony software, with the flight altitude set to 50 m. In addition, the overlap between the images was set to 70% and the movement speed of the UAV was set to 3 m/s. This software proved effective in mitigating the effects of the angle of solar radiation on the spectra by allowing the UAV to maintain a constant heading during the round trip. To ensure optimal conditions for data collection, all data were collected on sunny days between 9 a.m. and 12 p.m. local time to minimize atmospheric effects. Finally, the CSC hyperspectral data were smoothed using a Savitzky–Golay filter with a polynomial degree of 2 and a window size of 25 (4.2 nm on each side).

2.3. Methodology

2.3.1. Spectral Indices Extraction

To estimate the SM from optical data, vegetation indices are commonly used to reflect the condition of the plant cover at different SM contents. Since these indices can reflect the physiological state of plants, they serve as indicators of vegetation cover [50,51] and respond to changes in SM [52]. Therefore, several spectral indices, including the Normalized Difference Vegetation Index (NDVI) [53], the Normalized Difference Water Index (NDWI) [54], the Normalized Difference Moisture Index (NDMI) [55] and the Soil Adjusted Vegetation Index (SAVI) [56], were extracted from Landsat-8/9 and Sentinel-2 data using ArcMap (version 10.8.1). These indices were then analyzed for their correlation with in situ SM data.

2.3.2. Narrow Band Spectral Indices

Four indices, including the DI, RI, NDI, and PI were created using satellite bands and CSC bands. For this purpose, all possible two-band combinations were investigated in MATLAB R2023a using Equations (2)–(5) [31], and the correlations between all calculated indices and SM was evaluated.

$$D{I}_{\left(R_i,R_j\right)=R_i-R_j }$$

(2)

$$R{I}_{\left(R_i,R_j\right)=\frac{R_i}{R_j} }$$

(3)

$$ND{I}_{\left(R_i,R_j\right)=\frac{R_i-R_j}{R_i+R_j} }$$

(4)

$$P{I}_{\left(R_i,R_j\right)}=\frac{R_i-0.4401R_j-0.3308}{\sqrt{1+{0.4401}^2}}$$

(5)

where, R_i and R_j are the spectral reflectance of i and j bands, respectively.

2.3.3. Soil Moisture Modeling

In this study, the RF algorithm [57] was used for modeling SM. The RF algorithm is recognized as an ensemble learning method that is commonly used for regression and classification tasks. In the training phase, several decision trees are constructed. Each decision tree is trained on a random sample of the data with replacement. In addition, only a random subset of features is considered for splitting at each node of the decision tree, introducing randomness and diversity into the trees. During prediction, each tree in the RF independently generates its own prediction. For classification tasks, the class predicted by the majority of trees is chosen as the final prediction, while for regression tasks, the average prediction of all trees is taken.

The RF algorithm has been widely used in various soil studies related to soil properties such as soil texture and soil type, and its effectiveness in modeling SM with different datasets has been demonstrated in several studies [22,24,31,40,41,44]. The main advantages of RF in soil studies are its compatibility and high accuracy, especially for datasets with numerous predictors (features) and limited samples. The algorithm is characterized by its robustness to irrelevant features, and it eliminates the need for precise parameter settings [58]. Consequently, it is well suited for datasets with a large amount of data [59]. The ARTMO package [60] in MATLAB R2023a was used for SM modeling. After preprocessing the data, outliers and irrelevant data were excluded from the modeling process. The model included nine bands from Landsat-8/9, 11 bands from Sentinel-2, and 1252 bands from the CSC (covering the wavelength range between 420 and 850 nm).

In this study, the modeling was carried out independently for each dataset. Depending on the availability of datasets on the days of data collection, the total amount of in situ data for each dataset varied. Thus, a total of 213 data were collected for Sentinel-2, 324 data for Landsat-8/9, and 236 data for the CSC. Two-thirds of the in situ data were used to train the RF model for each dataset, and the remaining data were used to validate the trained models. Because of the small number of samples, a k-fold (k = 10) approach was used for cross-validation of the trained models. In this way, all the samples were used for training and validation. In addition, the results of the datasets on compatible days were compared with each other. Finally, the sensitive bands of each dataset to SM were identified, and uncertainty maps were created. Figure 4 shows the flowchart of the work.

2.3.4. Validation

In this study, four statistical indices, including the coefficient of determination (R²), the adjusted coefficient of determination (R²_adj), the root mean square error (RMSE) [61], and the mean absolute error (MAE) [62], were used based on the following equations.

$${\mathrm{R}}^2=1-\frac{{{\displaystyle \sum }}_{i=1}^n{\left(SM{C}_{ob,i}-SM{C}_{p,i}\right)}^2}{{{\displaystyle \sum }}_{i=1}^n{\left(SM{C}_{ob,i}-\overline{SMC}\right)}^2}$$

(6)

$${{\mathrm{R}}^2}_{adj}=1-\frac{n-1}{n-m-1}\left(1-{\mathrm{R}}^2\right)$$

(7)

$$\mathrm{RMSE}=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left(SM{C}_{p,i}-SM{C}_{ob,i}\right)}^2}{n}}$$

(8)

$$\mathrm{MAE}=\frac{{{\displaystyle \sum }}_{i=1}^n\left|SM{C}_{p,i}-SM{C}_{ob,i}\right|}{n}$$

(9)

where SMC_p and SMC_ob are the predicted and observed SM, respectively, $\overline{SMC}$ is the average of the observed SM, n is the number of soil samples, and m is the number of independent variables in the model that are used to predict the SM.

3. Results and Discussion

3.1. Effect of Soil Moisture on Multispectral Data

Table 1. Correlation of the bands and indices of the Sentinel-2 and Landsat-8/9 with soil moisture.

Features	Landsat-8/9			Sentinel-2
	Equation	R *	Band	Equation	R *	Band
Bands	y = −0.0144x + 33.4	−0.55	Aerosol	y = −0.007x + 81.01	−0.68	Aerosol
	y = −0.0113x + 31.7	−0.65	Blue	y = −0.0059x + 74.13	−0.68	Blue
	y = −0.0103x + 35.4	−0.65	Green	y = −0.0043x + 65.01	−0.68	Green
	y = −0.007x + 29.73	−0.64	Red	y = −0.003x + 52.69	−0.67	Red
	y = −0.0084x + 37.1	−0.61	Red Edge 1	y = −0.0009x + 23.98	−0.23	NIR
	y = −0.0036x + 23.9	−0.12	Red Edge 2	y = −0.0027x + 57.97	−0.67	SWIR-1
	y = 0.0018x + 0.71	0.08	Red Edge 3	y = −0.0026x + 49.72	−0.69	SWIR-2
	y = 0.0018x + 0.73	0.07	NIR	y = −0.0006x + 34.41	−0.30	TIR-1
	y = 0.0032x − 5.81	0.14	Red Edge 4
	y = −0.0101x + 51.3	−0.66	SWIR-1
	y = −0.0071x + 34.8	−0.64	SWIR-2
Indices	y = 54.345x + 6.56	0.57	NDMI	y = 65.453x + 7.94	0.43	NDMI
	y = 33.049x + 1.7	0.59	NDVI	y = 45.414x + 3.11	0.4	NDVI
	y = −42.71x − 2.94	−0.58	NDWI	y = −47.385x + 0.56	−0.32	NDWI
	y = 22.035x + 1.7	0.59	SAVI	y = 30.277x + 3.11	0.4	SAVI

* Pearson correlation.

shows the correlations between different bands and indices of Sentinel-2 and Landsat-8/9 satellites and SM. Of the Sentinel-2 bands, the SWIR-1 band showed the highest correlation with in situ SM data, with a correlation coefficient of 0.66 [23,63,64,65]. In addition, the Blue, Green, SWIR-2, and Red bands had correlations of −0.65, −0.65, −0.64, and −0.64, respectively, showing only minor differences in their rankings. Conversely, the Red Edge bands showed minimal correlation with SM, despite their indirect link to vegetation water content, which is sensitive to changes in SM [66,67]. While the Red Edge 3 band is pivotal for estimating crop water content, its correlation with SM was remarkably low in this study. This discrepancy can be attributed to the prevalence of bare soils and irregular irrigation schedules in the study area, and the finding differs from the results of Hegazi [23], who focused on a vegetated region. In addition, the time interval between irrigation and the response of plants to SM changes also influences these correlations [68]. In the case of the Sentinel-2 data, the SWIR-1 band showed a stronger correlation with the ground SM data than the SWIR-2 band, which is consistent with the results of Hegazi [23] and Hashim [69]. However, this contrasts with the results reported by Sanchez [70]. In contrast to the results of Zheng [71] and Hegazi [23], the Blue band in this study showed a remarkable correlation with the in situ SM data.

The results of the Landsat-8/9 data showed that of the Landsat bands, the SWIR-2 band had the strongest correlation with the in situ SM data (r = −0.69). The effectiveness of the SWIR band in accurately estimating SM has been emphasized in previous studies [71,72,73]. This is closely followed by the Blue band, with a correlation coefficient of −0.68, which contrasts with the results of Zheng [71] and Hegazi [23]. The SWIR-2 band showed a higher correlation with in situ SM data than the SWIR-1 band, which agrees with Sanchez [70] but differs from Hegazi [23] and Hashim [69]. The lowest correlation with the in situ SM data was found for TIR-2, with a correlation coefficient of 0.07. In general, the Landsat-8/9 thermal bands showed a modest linear correlation with SM, possibly due to their lower spatial resolution (originally 100 m, later converted to 30 m by processing). The considerable spatial variability of land use in the region, combined with the large pixel size, resulted in the presence of different land uses within a pixel, contributing to pixel non-uniformity and consequently reduced band performance. In addition, the different land uses in the region complicated the relationship between thermal bands and SM [74,75].

The correlation of both satellites with in situ SM data in the NIR bands was remarkably low due to water vapor absorption effects at a wavelength of 825 nm [76]. Moreover, the SWIR bands of both satellites showed similar performance [76] and are very suitable for estimating surface soil moisture in drylands [77]. For both satellites, the SWIR bands showed a higher correlation with the in situ SM data compared to the NIR bands [71,78,79].

As shown in Figure 5a, the NDVI values of the Landsat-8/9 were lower than those of the Sentinel-2 images for collected soil samples. Because of the presence of less green vegetation, the single band of Landsat-8/9 showed higher capability for SM retrieval, compared to Sentinel-2.

As shown in Table 1, the Sentinel-2 indices showed a higher correlation with SM compared to Landst-8/9. Apart from the NDWI, all indices showed a direct relationship with SM, meaning that as SM increased, the values of these indices also increased. The NDVI and SAVI were calculated using bands 4 and 5 of Landsat-8/9 and pair-bands 4 and 8 of Sentinel-2 images. The between-bands correlations for Landsat-8/9 and Sentinel-2 for the NDVI and SAVI were 0.28 and −0.38, respectively. The NDMI was calculated between pair-bands 5 and 6 of Landsat-8/9 and pair-bands 8 and 11 of Sentinel-2 images. The between-bands correlation for Landsat-8/9 and Sentinel-2 for the NDMI were 0.52 and −0.29, respectively. The NDWI was calculated between pair-bands 3 and 5 of Landsat-8/9 and pair-bands 3 and 8 of Sentinel-2 images. The between-bands correlation for Landsat-8/9 and Sentinel-2 for the NDWI were 0.4 and −0.32, respectively. The presence of more green vegetation in Sentinel-2 images can be proven by a positive Pearson correlation between pair-bands for Landsat-8/9-derived indices (Figure 5b) and a negative Pearson correlation between pair-bands for Sentinel-2-derived indices (Figure 5c). The between-bands correlation of Landsat-8/9 bands, achieved by increasing the data redundancy, reduced the performance of indices calculated from these bands for SM retrieval, whereas combining the pair-bands of Sentinel-2 led to improved accuracy of SM estimation.

3.2. Effect of Soil Moisture on Hyperspectral Data

In this study, bands with wavelengths below 420 nm and above 850 nm, which are characterized by noise, were excluded from the processing. The correlation between the different CSC bands and the in situ SM data is shown in Figure 6. The results showed that the effect of SM on different wavelengths is not the same. The lowest correlation was observed for the first bands in the 420 to 430 nm range. After this range, the correlation values increased sharply and peaked at the wavelength of 480 nm. Within the 480 to 770 nm range, the correlation of the different bands with the in situ SM data remained relatively constant, with the highest correlation coefficient recorded at 527 nm (r = −0.85). The correlation patterns observed in this study are consistent with the results of Jiang [80] and Mirzaei [21]. Beyond the wavelength of 770 nm, the correlation values gradually decreased.

3.3. Soil Moisture Estimation Using Indices

Figure 7 shows the 2D maps illustrating the coefficients of determination of the indices derived from the different bands of Sentinel-2, Landsat-8/9, and the CSC, for the estimation of SM. The results showed that the PI had an increased sensitivity to SM variations in all datasets [31,41]. Furthermore, the results showed that all CSC-derived indices performed better compared to the indices derived from Landsat-8/9 and Sentinel-2 data for SM estimation. In particular, the accuracy of the Landsat-derived indices exceeded that of the Sentinel-2-derived indices.

The optimal bands for each index, yielding the highest accuracy in SM estimation, are given in

Table 2. Optimal pair-bands for indices.

Index	Optimal Bands/Wavelengths			R2
	Sentinel-2	Landsat-8/9	CoSpectroCam (nm)	Sentinel-2	Landsat-8/9	CoSpectroCam
DI	Green, Red	SWIR-2 and TIR-2	427.5, 678.0	0.36	0.48	0.78
NDI	Red, Red Edge 1	SWIR-2 and TIR-2	525.4, 562.1	0.38	0.50	0.65
PI	Red Edge 1, SWIR-1	SWIR-2 and TIR-2	431.6, 523.6	0.42	0.48	0.79
RI	Red, Red Edge 1	SWIR-2 and TIR-2	525.4, 562.1	0.38	0.51	0.65

. In the case of the Sentinel-2 data, the PI created with Red Edge 1 and SWIR-1 had the highest accuracy in estimating SM, with an R² value of 0.42. This result confirms the claims of Liu [73] that the combination of Red Edge and SWIR bands significantly improves the accuracy of SM estimation. In the context of the DI, the optimal spectral bands for SM estimation were identified as Green and Red. In addition, the most effective band combination for creating the NDI and RI was identified as a combination of the Red and Red Edge 1 bands.

For the Landsat-8/9 data, the PI performed consistently better for SM estimation. In particular, the RI, which was formulated from SWIR-2 and TIR-2 bands, showed the best performance for SM estimation, with an R² of 0.51. The defined bands were also simultaneously identified as optimal for the DI, NDI, and PI and showed robust SM estimation abilities, with R² values of 0.48, 0.5, and 0.48, respectively.

In the CSC data, the PI showed higher sensitivity to SM changes compared to the other three indices. In particular, the PI achieved a remarkable R² of 0.79 in SM estimation when wavelengths of 431.6 and 523.6 nm were used. The optimal wavelengths for creating the DI, NDI, and RI were found to be 427.5 and 678, 431.6 and 523.6, and 525.4 and 562.1 nm, respectively. These wavelengths showed remarkable accuracy in estimating SM, yielding R² values of 0.78, 0.65, and 0.65, respectively.

3.4. Soil Moisture Retrieval Using Random Forest Algorithm

The data from all datasets were used as inputs to the RF model, and the results were validated with in situ SM data by 10-fold cross-validation (

Table 3. The accuracy of the random forest algorithm for soil moisture estimation using Sentinel-2, Landsat-8/9, and CoSpectroCam data.

Dataset	R	R2	R2adj	RMSE (%)	MAE
Sentinel-2	0.70	0.49	0.48	5.95	3.43
Landsat-8/9	0.81	0.66	0.65	4.61	2.42
CoSpectroCam	0.93	0.87	0.86	2.60	1.93

). The results underlined the superior performance of the CSC hyperspectral data in SM estimation compared to Sentinel-2 and Landsat-8/9 multispectral data [42]. In particular, the Landsat-8/9 data outperformed the Sentinel-2 data in SM estimation, partly due to the lack of a thermal band in Sentinel-2. The thermal bands in Landsat-8/9 had a high sensitivity to SM changes, contributing to its better performance. The R² for the model created with the Sentinel-2 data was 0.49, which, in contrast to the results of Liu [73], indicates a relatively modest potential for SM estimation; this discrepancy can be attributed to the limited size of the study area.

In contrast, the model built with Landsat-8/9 data yielded an R² of 0.66, confirming the significant potential of Landsat-8/9 data for SM estimation [26,81,82]. Impressively, the model created with CSC data achieved an R² of 0.87, highlighting the great potential of CSC hyperspectral data for SM estimation (See Supplementary Information, Figures S1 and S2 for SM maps predicted using Sentinel-2 and Landsat-8/9 data).

Figure 8 shows the predicted and observed SM levels for all datasets, with each date shown in a different color. The results showed a consistent trend wherein the models produced tended to overestimate the SM content in regions with less than 8% SM, while there was a tendency to underestimate the SM content in regions with more than 8% SM [41]. The degree of over- and underestimation of the SM was found to be more pronounced when using satellite data compared to when using the CSC. This can be attributed to the lower spatial resolution of satellites compared to the CSC, which leads to the non-uniformity and impurity of pixels [83]. A pixel may include both bare soil and vegetation while providing only a single value. Other factors contributing to the higher accuracy of CSC data are the higher spatial and spectral resolution and the use of preprocessing techniques such as noise removal and smoothing filters for CSC data.

SM maps derived from Sentinel-2 and Landsat-8/9 data on common overpass days (see Supplementary Information, Figure S3) revealed a consistent pattern in SM estimation for both satellites. Both satellites showed higher SM values in regions with vegetation compared to regions with bare soil, a phenomenon that can be attributed to the increased preservation of SM in vegetated areas.

Following the modeling of SM using satellite data, an uncertainty map for the SM estimation was generated based on the standard deviation of the observed and predicted values for each date (see Supplementary Information, Figures S4 and S5). The results indicated that the accuracy of both satellites decreased when estimating SM in vegetated areas, leading to an increased level of uncertainty in these regions. A comparison of the uncertainty map produced with Sentinel-2 and Landsat-8/9 data on days when both satellites had the same transit time over the study area (12 May 2022, and 21 June 2022) showed a higher degree of uncertainty in the maps derived from Sentinel-2 data. Nevertheless, the uncertainty patterns in SM estimation were similar for both satellites, with lower performance observed in areas with vegetation.

3.5. Band Sensitivity Analysis

The sensitive wavelengths of the CSC and the sensitive bands of the satellites in each iteration of the modeling process are shown in Figure 9. The findings showed that the Red Edge1 band of Sentinel-2 scored highest in 5 of the 10 replicates. In addition, the Red Edge 2 and Red Edge 4 bands secured the top position in three and two replicates, respectively. These results emphasize the importance and sensitivity of the Red Edge Sentinel-2 bands to variations in SM. However, as indicated in Table 1, the linear correlation between these bands and SM was comparatively lower. This discrepancy may be attributed to the inherent limitation of linear correlation in accurately representing the complicated relationship between these bands and SM, which was extracted non-linearly by the RF algorithm. The sensitivity of the Red Edge bands to SM changes is consistent with the findings of the study by Paloscia [73].

In the Landsat-8/9 data, the TIR-2 band secured the top position in nine of the ten replicates, highlighting its remarkable sensitivity to SM changes. However, as indicated in Table 1, the linear correlation between this band and SM was comparatively lower. The noticeable difference in sensitivity is probably due to the inadequacy of linear correlation in accurately describing the relationship between the dependent and independent variables. Since this relationship is not linear, it could be better elucidated by the RF algorithm. In the only remaining replicate, the Red band was identified as the most sensitive of the bands considered.

The analysis results of the sensitive bands in the CSC showed that the wavelengths most sensitive to SM variations are predominantly in the Blue and Red Edge ranges [31].

4. Conclusions

In this study, we aimed to fill the gap of previous studies on SM modeling by comparing the potential of multispectral and hyperspectral data for improved SM estimation. Our comparative analysis, using multispectral data from Sentinel-2, Landsat-8/9, and hyperspectral data from a CSC mounted on the Matrice 100 UAV, showed the superiority of hyperspectral UAV data for accurate SM modeling. This finding is significant as it enables the acquisition of high-resolution, site-specific SM information that is critical for optimizing irrigation practices and resource management in precision agriculture. This study not only contributes to a better understanding of SM distribution, but also highlights the potential of UAV-based hyperspectral data collection for improved SM monitoring. Despite the promising results in the accurate estimation of SM, there are still open questions that require further investigation. First, further insights into the influence of environmental factors such as land use change and climate variability on the accuracy of SM modeling need to be gained. Second, the integration of multiple datasets (e.g., the integration of UAV and satellite data) should be explored to improve the accuracy and spatial resolution of SM predictions. Third, investigating the potential of different machine learning algorithms for SM estimation could further improve modeling techniques. Finally, future research efforts should aim to evaluate the capability of radiative transfer models (e.g., MARMIT) for SM estimation.