Fusion of Reflected GPS Signals With Multispectral Imagery to Estimate Soil Moisture at Subfield Scale From Small UAS Platforms

This study proposes a low-cost and “proof-of-concept” methodology to obtain high spatial resolution soil moisture (SM) via processing reflected global positioning system (GPS) and a multispectral camera data acquired by small unmanned aircraft system (UAS) platforms. An SM estimation model is developed using a random forest (RF) machine-learning (ML) algorithm by combining features obtained from reflected GPS signals (collected by smartphones and commercial off-the-shelf receivers) in conjunction with ancillary vegetation indices from the multispectral camera data. The proposed ML algorithm uses in situ SM measurements acquired via SM probes as labels. A preliminary field experiment was conducted on 210 by 110 m (2.31 ha) crop fields (corn and cotton) in 2020 (from January to November, including crop planting through senescence time period) at Mississippi State University (MSU)’s the heavily instrumented North Farm to acquire data needed for the ML model to train and test. Our results showed that both fields could be covered by GPS reflectometry measurements with about 13 min of flight time at a 15-m altitude, and SM can be mapped with 5 × 5 m spatial resolution (corresponding to the elongated first Fresnel zone). The model is trained with and validated against eight in situ SM station datasets via tenfold and leave-one-probe-out cross-validation techniques. Overall, root-mean-square errors (RMSE) of 0.013 m <inline-formula><tex-math notation="LaTeX">$^{3}$</tex-math></inline-formula> m<inline-formula><tex-math notation="LaTeX">$^{-3}$</tex-math></inline-formula> volumetric SM and R-value of 0.95 [-] are obtained for tenfold cross validation. The proposed model reached an RMSE of 0.033 m<inline-formula><tex-math notation="LaTeX">$^{3}$</tex-math></inline-formula> m<inline-formula><tex-math notation="LaTeX">$^{-3}$</tex-math></inline-formula> and an R-value of 0.5 [-] in leave-one-probe-out cross validation. While having limited data, the results indicate that high-resolution SM measurement can be achieved with a low-cost GPS reflectometry system onboard a small UAS platform for use in precision agriculture applications.


I. INTRODUCTION
S OIL moisture (SM) content is one of the critical physical parameters for soil characterization [1], [2]. Accurate measurement of volumetric SM (VSM) is an essential part of site-specific precision agriculture (PA) management. Efficiently planned and managed irrigation water not only improves yields and product quality but also protects the environment and limited natural sources [3]. VSM traditionally can be measured using in situ SM probes in the field. This technique is generally regarded as accurate and reliable to the sensitivity of the SM probe and it is used as a reference for determining VSM [4]. However, high-resolution SM observations through SM probes can be time consuming, costly, and inefficient for large heterogeneous fields.
Over the last several decades, remote sensing approaches have become a popular way of retrieving spatially distributed VSM for regional and global studies. A variety of sensing techniques have been utilized to estimate VSM. Space-borne microwave remote sensing missions such as soil moisture and ocean salinity (SMOS) [5], soil moisture active passive (SMAP), Sentinel-1 [6], advanced SCATterometer [7], Fengyun-3B [8], and meteorological operational satellite [9] are typical examples of those that provide surface SM measurements on a large scale. These microwave-based satellite technologies mainly use L-band (1-2 GHz) or C-band (4)(5)(6)(7)(8) observations and produce a surface (top few centimeters) VSM product from active or passive observations. Their spatial temporal characteristics vary with respect to the sensing modalities. For instance, the European Space Agency's (ESA) SMOS and the National Aeronautics and Space Administration's (NASA) SMAP [10] are both dedicated SM missions that are instrumented with L-band passive radiometers and provide 36-km spatial and 1-3 days temporal resolution. Sentinel-1 is another ESA mission equipped with C-band synthetic aperture radar and can create SM maps with 1-km spatial and 6-12 day temporal resolution [6]. Although all existing satellite missions are critical for many large-scale climate, hydrometeorology, atmospheric research, and water resource management studies, they have a coarse spatio-temporal resolution that do not meet the practical requirements for many high-resolution subfield scale PA applications.
Unmanned aircraft systems (UASs) are widely available and utilized in PA, such as vegetation trait and stress monitoring [11], This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ [12], weed mapping [13], and irrigation management [14]. However, compared with these UAS applications, subfield scale VSM estimation studies are very limited and can be grouped as soil spectral reflectance [15], vegetation spectral reflectance [16], [17], and microwave response [18].
Soil spectral reflectance from UAS platforms can be, for instance, used for VSM estimation by characterizing its influence on spectral reflectance with laboratory measurements as studied by several previous works [19], [20]. These studies showed that reflectance level decreases as VSM increases due to a darkening of the soil surface color. Water Index SOIL [21], Shortwave Angle Slope Index [22], and Normalized Soil Moisture Index [23] are several spectral indices that are highly correlated with VSM. These indices have been tested in the laboratory at varying VSM levels over different bare soil types. The main drawbacks of soil spectral reflectance are that the performance depends on soil color and texture and organic material. In addition, the penetration depth is limited to the uppermost layer of soil. Furthermore, these indexes are mainly valid for bare or sparse vegetation areas and require atmospheric compensation for higher sensitivity. In a recent study [24], the authors showed the correlation between VSM and the brightness of UAS-based visible images for arid and semiarid regions. To minimize the effect of the atmosphere, the UAS flight height was kept close to the ground (1.4 m).
Establishing an inversion model of VSM based on vegetation spectral data is an alternative and well-studied approach for assessing SM. The principle is to use remotely sensed spectral observations strongly linked with VSM [25]. The studies show that many canopy spectral indices and biophysical and biochemical parameters are strongly correlated with VSM [26]. Recently some UAS applications have been performed using the vegetation reflectance technique [25], [27], [28].
The advantages of using the spectral index of vegetation can be summarized as sensitivity to target parameters and insensitivity to interference factors [29]. However, this technique depends on vegetation type and is not available for barren and sparsely vegetated areas. Additionally, changes in spectral indices can lag changes in climate/irrigation/precipitation. This lag may vary between 20 and 45 days depending on vegetation type [30], [31].
Using global navigation satellite system reflectometry (GNSS-R) for measurements of VSM is another attractive remote sensing technology. In this technology, all navigation satellites act as a transmitter, and a small size and low-cost receiver is sufficient to capture reflected signals from the ground. GNSS-R has been used in many geophysical applications, such as ocean surface winds [32], snowpack, sea ice [33], vegetation [34], and wetland [35] and flooding monitoring [36].
GNSS-R has already demonstrated the ability to sense SM from space-borne [37], airborne [38], and tower platforms [39]. For instance, many research groups developed models and algorithms to retrieve VSM from land observations of NASA's recent spaceborne GNSS-R mission called cyclone global navigation satellite system (CYGNSS) [40], [41], [42], [43], [44]. Most of these CYGNSS-based SM retrieval studies use peak power of delay-Doppler maps to calculate surface reflectivity. Although CYGNSS-based VSM estimations have high spatial and temporal resolution compared to other space-borne missions that measure VSM, its resolution under coherent scattering assumption is approximately 7 by 0.5 km with 1 ms integration time, and 25 by 25 km when noncoherent scattering is assumed [40]. For both cases, the achieved resolution by CYGNSS SM products is not enough for subfield scale PA applications. Airborne GNSS-R studies that create VSM estimations have been performed with special aircraft and GNSS systems. These may be expensive and produce spatial resolutions (30-300 m, depending on aircraft's speed and altitude) that are not fine enough for many PA applications [38], [45], [46].
To bring the GNSS-R technique into PA applications, we utilize a down-facing smartphone on a low-altitude (∼ 15 m) small UAS, where GPS signals reflected from the land surface are collected via the built-in GPS chip of the smartphone [47]. Due to the long pulse length of GPS signals (about 300 m for Coarse Acquisition C/A code), both direct and reflected signals interfere at low altitudes (e.g., below height of 100 m). Thus, we use a ground plane to block the direct signal and to isolate surface-reflected GPS signals that carry information around specular points (SPs). The preliminary results of the study indicate that reflected signals are highly correlated with surface features, including high reflectivity over ponds, sensitivity to land crop types, and moisture in the soil. Furthermore, the smartphone is capable of simultaneously acquiring multiple measurements from different surface locations (more than 10, depending on visible GNSS satellites during the experiments). Discussions on the UAS configuration and flight information, antenna and orientation effects, and the methods to prevent direct signals for GNSS satellites are provided in [47].
This article extends our aforementioned study over larger and more diverse datasets using both smartphone and custom GPS receivers along with auxiliary vegetation features from multispectral cameras to develop a data-driven high-resolution VSM retrieval algorithm. We train and test a machine learning (ML) algorithm using experimental data collected from UAS-based systems. The ML model aims to learn nonlinear and complex land surface characteristics between VSM and reflected GPS signal and multispectral data features. To this end, we conducted a preliminary field experiment on 2.31 ha crop fields (corn and cotton) from January to November in 2020 (including crop planting through senescence) at MSU's heavily-instrumented North Farm. The reflected GPS data are complemented by ancillary vegetation indices from a multispectral camera as well as in situ SM measurements.
The overarching objective of the present study is to develop a low-cost "proof-of-concept" solution for high-resolution (subfield scale) VSM retrieval from small UAS platforms. As a retrieval approach, we adapt our previously developed and validated space-borne GNSS-R VSM retrieval approach [48] to the geometries from small UAS platforms. To show the efficiency of the proposed model, we trained and tested the ML model independently using GNSS-R measurements for SM retrieval. Then, we trained the model with merged datasets. This study is fundamentally different from other UAS-based VSM retrieval systems by providing a GNSS-R microwave-based solution while the others focus mainly on image-based systems. The distinctions of this study can be summarized as follows. Fig. 1. Simplified structure of the proposed approach. UAS1a,b collects reflected GPS signals from the SPs via smartphone (for UAS1a) or U-blox GPS receiver (UAS1b). UAS2 collects multispectral images of the study field weekly to calculate vegetation spectral indices NDVI and VWC. Inputs C/N 0 , elevation, and azimuth angles from GPS receiver, and NDVI and VWC feed an RF model. The model trained with eight in situ SM measurements estimates SM for all study fields.
1) This study proposes a ubiquitous and low-cost microwavebased approach on a small UAS for subfield scale VSM estimations that could be suitable for PA. 2) This study fuses microwave measurements of reflected GPS signals with ancillary vegetation indices from a multispectral camera and develops a random forest (RF) ML model that can map VSM at high resolutions. 3) Contrary to the spectral reflectance-based technique (sensitive to greenness of vegetation), the proposed method efficiently operates for bare and vegetated areas. 4) Compared to space-borne reflected GPS data acquisition techniques, where direct and reflected signals are isolated in time, we use mesh metal planes to isolate reflected signal from the direct signal. This article is organized as follows. Section II describes the methodology, materials, data preparation process, and ML algorithm that were applied in the study. Section III presents and evaluates the experimental results. The implication of the findings and the potential future improvements of UAS-based VSM estimations are discussed in Section IV. Finally, Section V summarizes this study.

II. METHODOLOGY
This section describes the concept and overall structure of UAS-based VSM estimations from GNSS-R and ancillary measurements, material and hardware in developed systems and the flight details, the properties of the study field and details of the experiment, collected data types and preprocessing for the data, and details of the ML model for VSM estimation. Fig. 1 shows the concept and the simplified experimental steps for the proposed approach. A down-facing GPS antenna with a ground plane is used. The receiver can be a smartphone or commercial off-the-shelf receiver-UAS1a,b in the figure, and acquires reflected GPS signals over a specific predefined flight path (see Section II-C and Fig. 3) on each experiment date. This provides consistent data collection and near full coverage of the field. As UAS1 moves with a constant speed (e.g., 5 m/s in this work), the GPS receiver collects reflected carrier-tonoise density ratio (C/N 0 ) measurements from multiple SPs for each visible GPS satellite, simultaneously. Through these observations, we generate reflected C/N 0 maps of the field. Two different GPS receivers (i.e., built-in GPS chip and antenna of a smartphone, and a custom GPS receiver with linear polarized antenna) were used to acquire the reflected signals in the study. Experimental flights were performed for each receiver type within a 30 min lag.

A. Proof of Concept Overview
A separate UAS (UAS2 in the figure) is equipped with RedEdge-MX multispectral camera (https://micasense.com/ rededge-mx/), which provides spectral images in blue, green, red, red edge, and near-infrared (NIR) bands. The collected multispectral camera data are used to calculate specific vegetation indices, such as normalized difference vegetation index (NDVI), and vegetation water content (VWC), which are used as ancillary data in the developed ML-based VSM estimation architecture. UAS2, equipped with the multispectral camera, performs its weekly flight separately from UAS1. Hence, a weekly multispectral data are collected to observe changes in vegetation indices over the same study field.
To collect the ground truth SM measurements, eight HOBO SM probes are placed at different locations in the field. The placement of the probes is shown in Fig. 8. Each probe provides SM measurements for the top 5 cm VSM every 15 min. We acquire and record daily averaged surface SM values through the experimental period. These ground truth measurements are then used to train and validate the ML algorithm with the SP observations close to the probes. The ML model was trained using five features extracted from received GPS data and multispectral images. These features are reflected C/N 0 , elevation and azimuth angles from GPS data, NDVI, and VWC indexes calculated from multispectral orthoimages. The in situ SM values from the probes are used as the labels for the ML model. After the ML model is trained, the trained model is applied to test data, and an SM map of the study field is generated.

B. Materials and Hardware
In the beginning of the study, we started collecting reflected GPS signals with a smartphone. However, the preliminary tests showed that the irregular pattern of smartphone's bulit-in antenna made the interpretation for received signal challenging despite the use of a gimbal to fix its orientation. We decided to use an omnidirectional antenna with commercial GNSS receivers (see Section III-C and Table I) to minimize the impact of antenna and orientational variations on the measurements. Because of this, we started to fly a separate U-blox-equipped UAS on June 15th. We used a small-size quadcopter for U-blox-equipped UAS and a midsize quadcopter for smartphone-equipped UAS. The small-size UAS was built around a 550-mm carbon fiber frame (ZD550) with 13-in (33 cm) propellers and used Pixhawk arducopter as the flight controller, as displayed in Fig. 2(a). We used a U-blox application board C099 ZED F9P and AA.175 -Magma X linear polarized antenna as the GPS receiver. The antenna is attached underneath an 18-in (45.7 cm) circular mesh metal plate to block direct signals from GPS satellites affixed at the bottom of the UAS.
The midsize quadcopter was built on a Tarot 960-mm frame with 18 in propellers to carry smartphone-based (Xiaomi Mi8) GPS receiver. The receiver was attached underneath a 12-in (30.5 cm) circular mesh plate connected to a three-axis gimbal (Gremsy T3 V2) at the bottom of the drone depicted in Fig. 2(b). The built-in antenna of the smartphone is highly irregular, so the gimbal keeps the smartphone at a fixed orientation and level regardless of drone's pitch, yaw, and roll movement [49]. The metal mesh plate works as an electromagnetic shield that isolates the direct GPS signal from the reflected one. The UAS retracts its legs while collecting the reflected GPS signals. The U-blox GPS receiver is controlled with an onboard minicomputer on the UAS using a software (U-Center) designed for U-blox receivers while the smartphone uses an Android app named rinexON to record the signal. The GPS receivers on the UAS receive and record the reflected GPS signal while flying on the predefined flight path. A 5 in (12.7 cm) touch screen display is attached to the UAS for operating the custom GPS receiver. We also connected a downward-looking GoPro camera on the UAS for the real-time monitoring and recording of the field condition.

C. Study Field and Experiments
We conducted our study on a 210 by 110-m (2.31 ha) field (lower left corner 33 •   Foil Plant Science Research Center, Mississippi State University, Starkville, MS, USA during January-November in 2020. The mean annual temperature is 17.8 • C. July has the highest average monthly temperature, and January has the lowest. Annual precipitation is 1 329 mm. Average monthly humidity ranges from 55% to 66%. The elevation of the field is 92 m. At the time of the study, the field was organized with a split-plot arrangement and was planted with corn and cotton as main crops. The corn was planted on April 4 and harvested on September 11. The cotton was planted on May 4 and harvested on October 15. For all crops, the crop rows were spaced 38 in (96.5 cm).
The UAS follows a predefined flight path that contains nine north-south flight lines, as depicted in Fig. 3 Fig. 4(a) shows the flight dates performed with each UAS and measurements of eight SM probes (see Fig. 8 for the location of SM probes). Although we performed flights before corn planting and after the harvest of cotton, these flights were not shown in the figure. For the cornfield, SM probes' measurements started on May 7th, 2020 and ended on August 21st, 2020. SM probes in the cotton field started measurements on May 26th, 2020. Probe 2111 stopped its operation on September 15th, and probes 2113 and 2115 stopped on October 13th, 2020. After a sensor malfunction on August 20, SM data from probe 2125 were deemed unusable and removed. The measurements after August 2020 of this probe were removed from the dataset. Fig. 4(b) shows the statistics of each SM probe. The data corresponding to the in situ SM measurements period were used for the SM retrieval model development process.
The UAS flew at a constant speed of 5 m/s with one-second integration time. That means each 1 s observation corresponds to a trace of 5 m path on the ground. While the size of the first Fresnel zone at 15-m altitude for usable elevation angles (larger than 15 • ) is smaller than 5 m, 1-s observations were gridded into 5 × 5-m pixels and a VSM estimation is done for each resolution grid, where GPS observations are acquired.

D. Data Preparation
The reflected GPS (L-band) signal power is not only a function of VSM but also surface roughness and vegetation. Vegetation cover affects reflected signal power. Vegetation cover on a soil layer attenuates incident and reflected signals and adds its reflection [50]. The effect of the vegetation cover should be considered during the model development process. The vegetation cover can be characterized using several different indices. NDVI is an often-used index that shows vegetation activity [51]. VWC is another essential parameter used in retrieving VSM via microwave remote sensing. The weekly NDVI is calculated from light reflectance in the red and near-infrared wavelength bands in multispectral images as VWC is calculated by using NDVI and the stem factor s f of the vegetation [52] as where NDVI max and NDVI min are the annual maximum and minimum NDVI at a given location, respectively. The stem factor is an estimate of the peak amount of water residing in the stems, and for crops, it is taken as 3.5 by following SMAP's retrievals approach [53]. NDVI and WVC maps at 8-cm spatial resolution are aggregated to 5-m by averaging to match the GNSS-R resolution. We have used the closest available multispectral data to add NDVI and VWC to the GNSS-R features. The recorded GPS navigation messages are decoded, and C/N 0 ratio, elevation, and azimuth angles are extracted for each visible GPS satellite during the flight. The C/N 0 is the ratio of received carrier power to noise density in the receiver system. The elevation describes the angle of a satellite relative to the horizontal plane. The azimuth is the angle between the satellite and true north. The observations that correspond to elevation angles lower than 15 • were masked since the signal strength is much smaller and fluctuates more due to possible multipath as a blend of direct and reflected signals is expected at low elevation angles. In addition, the C/N 0 observations smaller than 15 dB were removed from the dataset. Specular reflection points are determined using simple image theory for each reflected GPS signal using the position of the UAS and the GPS satellites. During the study, only SP reflection is assumed and surface roughness effects are ignored. Each specular reflection point and its corresponding data are gridded into 5 × 5 m to generate a daily input map.

E. Machine Learning-Based VSM Estimation and Validation
This section describes the development of ML-based model for VSM estimation, as illustrated in Fig. 1, utilizing the collected dataset described in Section II-D. In this study, an RF regression algorithm was utilized due to its higher performance in comparable satellite-based VSM estimation studies [43], [48], [54], [55]. The RF model is trained and validated with a total of five input features C/N 0 , elevation angle, azimuth angle, NDVI, and VWC. in situ SM observations are used as the reference SM labels for the RF model. The required ground truth VSM labels were provided by eight SM probes located in the field acquiring data throughout the experimental period.
The developed RF model in this study contains ten trees with a maximum split size of six for each tree. The RF model is trained using a least-squares boosting ensemble strategy with a learning rate of 0.75. During the model training, we used only samples whose SP location was close (< 5 m) to any SM probes. Two different training and cross-validation procedures are used to assess the performance of the proposed SM model. These are tenfold cross validation and leave-one-probe-out validation. For the tenfold cross validation, the dataset is first split into ten folds, then nine folds are used as the training set, and the remaining fold is used as the testing set. The overall evaluation result is the averaged result of each fold. In the second cross validation, we have used the leave-one-probe-out approach. This leave-one-out cross validation is one of the effective ways to evaluate how well a model will perform on out-of-sample data.
In the leave-one-probe-out approach, data for a single SM probe are used as the test dataset while data for all other SM probes are used for training. Having a dataset of eight SM probes, seven were chosen for the training model, and the remaining probe's data were used for the test. These procedures were repeated for each SM probe data. The developed ML model performances are evaluated using root-mean-square error, unbiased RMSE, and R-value performance metrics.

III. RESULTS
The results will be given in two parts: first, for a selected experiment day, the maps that belong to measured features and the estimated VSM map will be provided. Second, the cross-validated performances will be summarized in the terms of performance metrics.

A. Experimental Results
A total of 50 flights were performed with smartphoneequipped UAS during January 7 to October 15. During these flights, the smartphone collected data from an average of 10.14 (standard deviation ±1.35) GPS satellites. After masking lowreflected power signals and low elevation angle satellites, data from an average of 6.28 (±1.05) GPS satellites remained for VSM estimation each day. From the utilized GPS satellites, an average of 4736 (±878) SP reflections were recorded daily. For the U-blox GPS receiver-equipped UAS, 40 flights were performed between June 15 and November 12. With these measurements, an average of 7.56 (±0.8) usable GPS satellites and 5378 (±629) SP observations were recorded. The results show that the SP measurements with smartphone receiver cover 84% of the study field, while the experiments with the U-blox GPS receiver cover 86% of the field on average. It is also possible that multiple SPs fall into the same 5 × 5-m grid cells. We observed that 71% and 82% of the cells have multiple specular measurements for smartphone and U-blox GPS receivers, respectively. Table I reports the statistics about the performed experiments.  VSM estimations are obtained from the developed ML model for each GPS observation. These individual VSM estimations for the same grid were averaged and assigned as the final VSM estimate of the grid for that day. The results show that the current flight duration of 13 min is sufficient to acquire necessary measurements to cover more than 80% of the entire study field. Fig. 7(a) and (b) shows the estimated VSM using U-blox GPS receiver and smartphone on July 10, and Fig. 7 (c) shows the differences between two VSM estimations. Although the proposed method provides VSM estimation for most parts of the field, there are grids without SM estimation. Uncovered parts generally are the east and west ends of the field, and they are parallel to the flight path. This result indicates that the flight lines should start before the boundary of the study field. From the figure, it is clearly seen that there are a smooth changes over the field and the overall VSM for corn (upper part) and cotton (lower part) are slightly different. This difference and the alley between the two crop fields are more noticeable in the U-blox VSM estimation map. The figure shows that the overall VSM values via smartphone are slightly higher than VSM estimations via the U-blox unit.

B. Performance Evaluation of the VSM ML Model
In this section, we evaluated UAS-based VSM predictions against the in situ VSM measurements from eight SM probes on the field. Fig. 8 shows the positions of SM probe stations. The evaluation was performed with two different training-validation strategies: tenfold cross validation and leave-one-probe-out validation. Fig. 9 Fig. 10 shows the estimated VSM for leave-one-probe-out cross-validation and in situ measurements in time series. This validation process is more challenging because the ML model that is tested over measurements at a single probe has been trained with measurements from different portions of the field without any observations from the test probe. In this challenging approach, overall RMSE and R-value were calculated as 0.033 m 3 m −3 and 0.50, respectively. ubRMSE was calculated between 0.016 and 0.037 m 3 m −3 for different in situ stations.
Moreover, we independently trained and tested the SM retrieval model using GNSS-R features. The results indicated that the fusion of spectral and GNSS-R features provides better SM retrievals for both devices in both cross validations. Using only GNSS-R features, the smartphone dataset reached an R-value of 0.38, whereas the U-blox dataset reached 0.48. By adding spectral features, NDVI and VWC, to the dataset, we have reached R-values of 0.48 and 0.50 for smartphone and U-blox, respectively. The spectral inputs provided information related to the vegetation to correct the effect of vegetation in the proposed method. Table II summarizes the overall performance of the proposed model for different cross-validation techniques and independent comparison of spectral and GNSS-R features.

C. Comparison of Smartphone and U-Blox GPS Receivers
The overall performance results of the study show that a UAS-based smartphone or a specially designed GPS receiver could be used to collect reflected GPS signals. Highly accurate VSM estimations can be achieved from data collected by both receivers with the proposed ML-based approach. However, as seen in Table II, the U-blox GPS receiver data model has comparably better performance than the smartphone receiver data model. The receive antenna for the U-blox GPS receiver is more uniform than the smartphone antenna.
To make a better comparison, we performed two UAS flights with 15-min lags for most of the study period. First, we flew the smartphone-based receiver, and then, we flew the U-blox-based receiver. The satellite's elevation angle did not change significantly between these two flights as it takes only 13 min to complete the full flight path. The received C/N 0 for both U-blox and smartphone receivers and the elevation angle during both flights are shown in Fig. 11. The red dotted lines represent the reflected signal received by the smartphone in dB-Hz, while the black dashed line represents the signal received by the U-blox GPS receiver in dB-Hz. The solid blue line represents the satellite's elevation angle in degree. Due to the near constant elevation   Fig. 3. While both receivers acquired similar and strong reflected signal over the water body, it can be seen that the U-blox receiver exhibits higher variation (and potentially higher sensitivity) over the field compared to the signal strength acquired by the smartphone. In addition, the results in Table I show that smartphone collected data from fewer GPS satellites than the U-blox GPS receiver. These observations are confirmed with the developed ML model performance results provided in Table II, where the ML model developed for smartphone is comparably less accurate than the U-blox model. This indicates the higher reliability of U-blox receiver over the smartphone. This could be due to the GPS antenna pattern of the smartphone being more complex as the received signal power changes along the antenna pattern at different elevation angles. Although we have used a gimbal to obtain a consistent smartphone orientation during the experiments to reduce this effect and incorporated elevation angles into the ML model, the effect of orientation and antenna issues was observed in the final results.

IV. DISCUSSION
This section will discuss the relationship between the experimental flight plan, collected data, developed ML model, and resultant VSM estimations. In the study, a 13-min flight duration with a 15-m flight altitude allowed the sensor to collect enough data to cover more than 80% of the study field (about 2.31 ha) with 5 × 5 m grid sizes. Decreasing the flight speed and/or altitude will allow for more data at finer spatial resolutions at the expense of longer flight times. However, longer flight times require higher battery capacity, which adds more complexity to the UAS design. Also, decreasing flight speed leaves the flight more susceptible to wind and motion, which, in turn, affects the quality of observations. The selected flight parameters were chosen to balance flight duration, stability, and field coverage.
Averaging multiple reflection data points over a single grid (in this study it was 5 × 5 m) can improve performance of VSM estimations. However, reducing the speed of UAS makes the stability of the UAS more difficult and will change the GPS antenna orientation that results in changes in measured reflection power. As an alternative, the UAS can repeat the flight path and more reliable data can be obtained. Although, in this study, we recorded the orientation (roll, pitch, and yaw angles) of the UAS and the GPS antenna, these parameters are not used in development process of the ML model due to the lack of ground truth data to prevent possible overfitting of the model.
The independent evaluation of GNSS-R shows the importance of GNSS-R features for SM retrieval. Merging GNSS-R features with spectral features provides an improvement in SM estimation. Although the results of the fusion of GNSS-R and spectral features have promising improvements in SM retrieval, it is important to note that we are performing this analysis with a limited number of data. A more detailed analysis will be  performed in future studies covering more frequent data collection. We want to note an essential point for GNSS-R-based SM retrieval; the reflected signal (C/N 0 ) is not only a function of SM but also antenna pattern, GPS satellite output power, elevation angle, surface roughness, and attenuation effects of vegetation. Although the output power of GPS satellites is consistent in time, this power level is not the same for each GPS satellite.
The distance between the GPS satellite and the study field (a function of the elevation angle) is another critical parameter that affects measured C/N 0 . Normalizing measured C/N 0 with GPS satellites' transmit power and elevation angle may help increase the SM retrieval accuracy. In this study, we used GPS data only. For future research, we will consider other GNSS missions, such as Galileo, GLONASS, and Beidou for the SM retrieval. Although surface roughness affects on L-band reflection, it is neglected in this study. If surface roughness parameter can be measured via another drone equipped with a LIDAR, it may be included in the ML model for enhanced VSM estimation. We filtered the reflections with elevation angles lower than 15 • since the observation may belong to incoherent reflection. This threshold level could vary depending on the surface conditions. Increasing this threshold level will reduce possible incoherent reflections. However, remaining GPS reflections may not be enough for complete spatial coverage of the study field.
In this study, we have used NDVI and VWC to measure vegetation information. However, other vegetation indices, such as the green NDVI and vegetation health index, can be considered. To obtain vegetation indices, a separate UAS with a multispectral camera is utilized. Both the GPS receiver and the camera can also be equipped in the same UAS. This would allow a single UAS to generate a VSM map of the field. Current experiments are done with a separate UAS to obtain more stable data, but our future work will include the development of a joint GPS receiver and multispectral camera load for a single UAS platform. As another alternative, instead of a multispectral camera, a low-cost RGB camera can be used, and different vegetation index can be calculated as input to the VSM retrieval model. Detailed analysis of such possible changes should be done in order to determine the optimal combination of sensing strategy for UAS-based high-resolution VSM mapping.
Finally, this is an initial study that was performed with limited flight and in situ measurements. In particular, the limited number of SM probes in the field makes the developed ML model's spatial evaluation difficult. As seen in Fig. 8, some SM probes are located very close to each other, so the spatial variation of reference SM is limited. In addition, temporal SM variation of the field is limited between 0.015 and 0.035 m 3 m −3 . Although there were multiple rain events during the study period, we could not collect data on rain days because of difficult flight conditions for the UAS. In addition, we used daily averaged SM probe measurements from the field as reference SM. However, this makes it challenging to catch diurnal changes and negatively affects SM retrieval performance. To remedy these drawbacks, one can use a mobile robotic instrument that can traverse crops to take SM measurements over an entire field coincident with UAS flights. This can allow us to train and test the model spatially at more locations. In addition, we plan to increase the number of SM probes in the study field for a better site-independent model. Also, we plan to conduct more flights over longer periods that cover at least two growing seasons to increase the time independence of the model.

V. SUMMARY AND CONCLUSION
This article illustrated the UAS application of a low-cost and high-resolution (subfield scale) GNSS-R-based remote sensing methodology for surface VSM in PA. Surface-reflected GPS signals, in combination with vegetation indices via a multispectral camera, provide sufficient modeling information to estimate surface SM accurately. In the study, nonlinear and complex land surface characteristics between VSM and reflected GPS signal and multispectral data were modeled using an RF ML algorithm that uses limited features, C/N 0 , elevation and azimuth angle, NDVI, and VWC. VSM estimation was performed with a developed ML model for the crops' complete growing season. Surface VSM estimation was accomplished with an overall RMSE of 0.013 m 3 m −3 and correlation of 0.95 for tenfold cross validation. The obtained performance score of the model using different cross-validation techniques shows the capability of the application for different crop fields. The proposed method can enable accurate prediction in the agricultural application to identify dry/wet spots and water-stressed crops for different time periods. A more accurate and better spatial resolution could be obtained by adjusting different flight plans (flight altitude, flight speed, and flight duration). We consider this work as a preliminary "proof of concept" study and future research will use more in situ SM measurements and flights to obtain an enhanced-trained ML model that can allow more detailed performance analysis of the proposed approach.