A soil moisture‐based framework for guiding the number and location of soil moisture sensors in agricultural fields

Soil moisture information is a key variable for guiding in‐season management decisions in rainfed and irrigated agricultural systems. However, methods for deciding the number and location of soil moisture sensors (SMS) per field still remain poorly explored in the scientific literature. The goal of this study was to evaluate a quantitative framework based on soil moisture‐based management zones (MZs) to determine the minimum number and tentative deployment location of SMS. Multiple spatially intensive (n > 100 observations) surveys of near‐surface (0–12 cm) soil moisture were conducted during the fallow periods and early growing seasons of 2017, 2018, and 2019 on three agricultural fields using a calibrated handheld soil water reflectometer. The fuzzy C‐means (FCM) clustering method was used to delineate MZs based on the soil moisture surveys, and the silhouette clustering evaluation method was used to identify the optimal number of MZs per field. Then, a sensor location index that considered the distance to the MZ boundaries and the FCM membership grade was developed to identify the tentative optimal deployment location of SMS. The proposed method effectively identified field areas with distinct soil moisture regimes and revealed the complex soil moisture spatial patterns that were not captured with elevation or soil texture alone. Dividing the fields using soil moisture‐based MZ reduced the intrazone soil moisture spatial variability by about 50% compared with that of the entire field. In the three studied fields, a total of two SMS were sufficient to capture the salient soil moisture spatial regimes.


INTRODUCTION
In situ soil moisture sensors (SMS) provide field-specific and timely information that can be used to guide in-season management decisions in rainfed (e.g., amount of fertilizer and seeding rate) and irrigated (e.g., irrigation scheduling) agricultural systems. In the United States, the adoption of soil moisture sensing technologies remains low, with only ∼11% of the U.S. farms using soil moisture sensing devices (Kukal et al., 2020). However, with the increasing need for sustainable use of water resources and the advent of more affordable and reliable soil moisture sensing technologies, the adoption of SMS by producers, irrigation managers, and scientists is expected to increase in the near future. New SMS offer innovative sensor designs, practical installation methods, and convenient on-board wireless data telemetry (e.g., radio frequency, cellular networks, and local long-range wide-area Vadose Zone Journal networks) that adapt to a wide range of consumer needs and that enable the deployment of local field-level soil moisture networks for precision irrigation scheduling (Hedley et al., 2012). However, a challenge that remains insufficiently addressed in the scientific literature resides in identifying the optimal number and position of a limited number of soil moisture monitoring locations in agricultural fields, a problem that largely depends on discerning field areas with clearly defined and temporally stable soil moisture spatial patterns. A major consideration when deploying an array of SMS is the limited scalability of the resulting monitoring system. Installing SMS in agricultural fields often requires manual soil trenching or soil augering, the additional deployment of hardware for data collection, and data telemetry that can conflict with frequent farming operations, and often involves additional soil sampling to develop a soil-specific sensor calibration. As a result, the installation of SMS is often approached as a semipermanent setting (i.e., several growing seasons), where delineation of zones with more homogeneous soil moisture conditions relative to that of the entire field is a cost-effective alternative. Traditionally, the delineation of field management zones (MZs) has been aimed at guiding the application of fertilizers based on spatial patterns of grain yield and surface soil chemical properties (i.e., soil pH, P, organic matter; Schepers et al., 2004), apparent soil electrical conductivity in combination with soil texture and elevation (Fraisse et al., 2001;Peralta et al., 2015;Reyes et al., 2019), yield maps collected over multiple growing seasons and different phases of the crop rotation (Basso et al., 2007), and using proximal and remote canopy reflectance indices (e.g., normalized difference vegetation index [NDVI]) from unmanned aerial vehicles (Corti et al., 2020;Ohana-Levi et al., 2019) and satellite platforms (Boydell & McBratney, 2002). In recent years, advancements in variable rate irrigation technology coupled with the increasing need to optimize irrigation efficiency to conserve surface and groundwater resources has increased the focus on delineation of MZs with emphasis on soil moisture (Haghverdi et al., 2015;Hedley et al., 2012;Sadler et al., 2005). However, existing research in this area heavily relies on proxy variables for soil moisture rather than actual soil moisture observations. For instance, delineation of site-specific MZs for irrigation scheduling typically includes clustering techniques based on multiple combinations of apparent soil electrical conductivity, elevation, topographic wetness indices, and soil texture (Boluwade et al., 2015;Hedley & Yule, 2009;Moral et al., 2010;Reyes et al., 2019); remote sensing vegetation indices (e.g., NDVI; Haghverdi et al., 2015;Reyes et al., 2019); and characterization of soil water-holding capacity (Hedley & Yule, 2009;Oldoni & Bassoi, 2016;Zhao et al., 2018). Previous studies at the watershed level have shown that properties dominating the soil moisture spatial variability change with soil moisture conditions. For instance, during wet

Core Ideas
• Field surveys of near-surface soil moisture were used to delineate management zones. • Soil moisture-based management zones were used to find the number and location of sensors. • Two soil moisture-based management zones captured the salient soil moisture spatial patterns. • In situ soil moisture observations revealed spatial patterns not evident with proxy variables.
conditions soil moisture patterns are mainly driven by soil porosity and hydraulic conductivity, but during dry conditions the major controlling factors are topography, slope aspect, and clay content (Famiglietti et al., 1998). Thus, because in situ soil moisture observations reflect the combined interplay between extrinsic landscape properties, management practices, and intrinsic soil physical properties, delineating MZs using in situ observations has the potential to reduce the uncertainty associated with proxy variables in applications where soil moisture information is needed. Therefore, we hypothesize that the delineation of MZs directly using soil moisture observations will elucidate field spatial patterns that cannot be detected using proxy variables. The goal of this study was to test a quantitative soil moisture-informed clustering framework to determine (a) the minimum number of SMS per field and (b) the tentative deployment location of SMS within each agricultural field. The framework is presented as a case study for three agricultural fields in which we conducted a series of spatially intensive soil moisture surveys to identify the salient soil moisture patterns.

Experimental fields
The study was conducted in three agricultural fields located in central Kansas during 2017 (Field A), 2018 (Field B), and 2019 (Field C) (Table 1). Most surveys were conducted during the fallow period, and a few surveys were conducted in fields that were planted with crops in early growth stages. Surveys in fields with planted crops were only conducted within 3 wk of planting when crops had little influence in the soil water balance ( a Letters were assigned to fields to protect the identity and location of the producers. b Sampling grid size is the approximate average size for all surveys, which was defined based on field area, available workforce, and a minimum number of samples required to interpolate soil moisture using a 5-m grid resolution. of deep, well-drained, alluvial soils with slopes <1% corresponding to the silt clay loam Detroit series (fine, smectitic, mesic pachic Argiustolls) and silt loam Hord series (fine-silty, mixed, superactive, mesic cumulic Haplustolls) (Soil Survey Staff, 2020). Field B is located near the city of Hutchinson, KS, and has an area of 58 ha. The cropping system of Field B consists of a minimum tillage operation based on irrigated corn and soybeans under a center pivot. This field is characterized by a rolling terrain dominated by deep loamy alluvium soils corresponding to the Avans series (fine-silty, mixed, superactive, mesic udic Argiustolls) with 2-4% slope and aeolian sandy loam soils corresponding to the Saltcreek series (fine-loamy, mixed, superactive, mesic udic Argiustolls) with 3-5% slope (Soil Survey Staff, 2020). The rolling terrain has a strong impact on the velocity and direction of runoff water flow in this field, which leads to areas of the field with pronounced rill erosion and localized accumulation of stubble after copious rainfall events. Field C is also located in the central portion of Kansas near the town of Moundridge, KS. The field has an area of 22 ha, and the cropping system is based on continuous no-till irrigated corn using a center pivot. This field has an upland area characterized by deep and moderately well-drained silt loam soils with <1% slopes that belong to the Crete series (fine, smectitic, mesic pachic udertic Argiustolls). The bottom slope is characterized by the Farnum series (fine-loamy, mixed, superactive, mesic pachic Argiustolls) with loam soils and slopes ranging from 1 to 3% and with the presence of fine gravel in the top soil horizon (Soil Survey Staff, 2020). The region encompassing the three fields has an approximate average annual rainfall of 800 mm, a mean annual temperature of 13˚C, and mean minimum and maximum air temperatures of −1 and 27˚C, respectively (Kansas Mesonet, 2021;Patrignani et al., 2020). The region belongs to the Köppen climate class with humid continental hot summers with year-round precipitation.

Soil moisture surveys
Soil moisture was measured in the top 12 cm of the soil profile using a handheld soil water reflectometer (HydroSense II CS659, Campbell Scientific). The Hydrosense is a lightweight portable sensor that consists of two 12-cm-long stainless steel rods attached to an epoxy sensor head. The sensor is connected to a display equipped with an onboard GPS receiver (±3 m accuracy) that is used for logging the soil moisture readings, the timestamp, and the associated geographic coordinates of each observation. The sensor also stores raw variables such as the apparent dielectric permittivity, which was used to develop a custom sensor calibration equation for accurate determination of volumetric water content. Field measurements of soil moisture were concentrated during the fallow periods and early crop stages to minimize the interference of frequent irrigation events and actively growing vegetation in the determination of soil moisture patterns (Hupet & Vanclooster, 2002). Fallow and early growing season periods also facilitated the manual collection of soil moisture observations during the intensive surveys, a task that would have been more labor intensive and time consuming in the presence of standing crops and frequent irrigation events. The timing of the soil moisture surveys was dictated by weather and field soil moisture conditions. For instance, the surveys on 5 and 8 Apr. 2019 were conducted before and after a rainfall event, respectively. The range of soil moisture conditions that we were able to capture was somewhat limited by the weather conditions at each location during the fallow periods. Each intensive soil moisture survey consisted of collecting georeferenced observations in the top 12 cm of the soil profile following a rough grid (Table 1). In total we conducted 16 soil moisture surveys across the three agricultural fields from May 2017 to June 2019, totaling 4,659 observations (Table 1, Supplemental Figure S2). A custom sensor calibration was developed in laboratory conditions using bulk soil from different areas of the fields. The bulk soil was first air dried, then ground to pass a 2-mm sieve, and then packed into cylindrical containers (2,980 cm 3 ) to create soil columns with known volumetric water content. A linear model (Ledieu et al., 1986) was developed relating the apparent dielectric permittivity measured by the sensor with the observed volumetric water content of the containers: T A B L E 2 Statistical properties of the soil moisture surveys collected during the fallow (F) periods and early growing season (EG) at the three agricultural fields including survey mean, CV, minimum and maximum soil moisture observations, best fitted semivariogram (SV), semivariogram range, the nugget-to-sill ratio (NSR) of the selected semivariogram, and the cumulative antecedent precipitation for each survey  The survey on 3 May 2019 only had 69 observations and the semivariogram was ill-conditioned, with a range value far exceeding the size of the field. A minimum of 100 points are typically required for stable kriging interpolation (Saito & Goovaerts, 2000).
where θ (m 3 m −3 ) is the volumetric water content, a (unitless) is the real part of the apparent dielectric permittivity, and the two empirical coefficients were determined by fitting a linear model. The custom calibration equation resulted in a RMSE of 0.027 cm 3 cm −3 and r 2 = .97 on the calibration dataset with volumetric water contents ranging from 0.013 to 0.448 cm 3 cm −3 (Supplemental Figure S1). For reference, the factory default calibration resulted in a RMSE = 0.041 cm 3 cm −3 on the calibration dataset. Field-level soil moisture for each intensive soil moisture survey was estimated using the ordinary kriging spatial interpolation method. The first step of the spatial interpolation process consisted of detrending the field observations collected in each specific survey using a quadratic linear regression model based on field geographic coordinates. The second step consisted of computing the regression residuals by subtracting the regression model from the soil moisture of each point in the survey to generate an empirical semivariogram. Spherical, exponential, and Gaussian semivariogram models were fitted using ordinary least squares. The semivariogram model with the lowest RMSE was selected to represent the spatial dependence in the interpolation step. The exponential semivariogram model resulted in the best fit in 63% of the surveys (10 out of 16 surveys), whereas the spherical semivariogram model resulted best in the remaining 37% of the surveys (6 out of 16 surveys) ( Table 2). The fitted semivariogram for each survey was used to interpolate the detrended field soil moisture observations using ordinary kriging, and then the spatial trend was added back to the interpolated detrended soil moisture. For the spatial interpolation, we discretized the area of the field using a regular grid of 5-m (i.e., grid cell area of 25 m 2 ) individual grid cells. The grid spatial resolution was defined so that the fifth percentile of the distances between sampling points of a given survey was represented by at least four grid cells (Hengl, 2006). The kriging routine was implemented in Matlab and was used with a local search neighborhood with a maximum of 500 m and a maximum of 15 neighboring observations. Since in this study soil moisture observations between different surveys were made at different positions within the field, differences in the soil moisture spatial patterns between different surveys contain an additional source of measurement variation. To eliminate this measurement variance, observations would need to be made at the same exact location across surveys, which can be impractical in large production fields and in surveys with large number (>200) of observations.

Determination of the number of MZs
A common approach to delineate MZs based on spatial features for agricultural applications is the use of unsupervised clustering. Fuzzy C-means (FCM) (Bezdek, 1981) is a widely used clustering technique to partition heterogeneous agricultural fields into more homogenous MZs (Peralta et al., 2015;Reyes et al., 2019;Yari et al., 2017). In contrast with other unsupervised clustering techniques like K-means (MacQueen, 1967), a distinct advantage of the FCM method is that any grid cell of the grid representing the area of the field can be assigned partial membership to different clusters (i.e., soft clustering) instead of assigning grid cells to mutually exclusive clusters (i.e., hard clustering). In other words, each discrete grid cell representing the field is assigned a membership grade in the continuous range from zero to one according to the degree of similarity with other grid cells of the same cluster (in this study a cluster represents a MZ). A value of one represents perfect membership and a value of zero represents complete lack of membership to a specific cluster. This unique feature typically makes FCM the preferred method to characterize continuous variables like soil moisture (Reyes et al., 2019), soil physical and chemical properties (Peralta et al., 2015), and crop yield (Burrough, 1989). The continuous range of the membership grade also enables integration of the continuous membership grade into a framework with additional continuous covariates. The FCM clustering was implemented using the Matlab function ffcmw from the Matlab File Exchange (Cococcioni, 2020). Because the FCM clustering method requires specifying the number of clusters as input, we combined the FCM clustering method with the silhouette method (Rousseeuw, 1987) to objectively evaluate the optimal number of MZs per field given a set of normalized variables. The silhouette method computes the mean separation Euclidean distance between all grid cells classified within a cluster and evaluates the interclass similarity among grid cells in other clusters, hence the resulting silhouette value serves as an objective metric to define the optimal number of zones per field. Silhouette values range between −1 and 1, where the highest value indicates the optimal number of clusters. In our case, silhouette values were computed for a range between one and six clusters, which spans the typical range of site-specific MZs in agricultural fields. In the special case in which the end goal is not to deploy the least number of SMS, but rather deploying a predefined number of SMS (e.g., following budget constraints), the step involving the silhouette method can be skipped and a predefined number of clusters can be enforced in the FCM clustering technique.
Before using the interpolated soil moisture surveys as input variables in the FCM, interpolated grids of soil moisture for each survey were normalized relative to the field average, following the same approach commonly used in time stability analysis (Vachaud et al., 1985): where β is the soil moisture relative difference at location i at time j, (θ ) is the expected volumetric water content in location i at time j obtained from the spatial interpolation (i.e., the resulting raster map of soil moisture for a given survey), and θ is the average field volumetric water content at time j.
To contrast the classification of MZs between in situ soil moisture observations and proxy variables, we also included in the analysis auxiliary variables such as the percentage of sand content, percentage of clay content, and elevation that were also normalized using Equation 2 before the clustering analysis. The similarity of the delineated MZs solely based on in situ soil moisture surveys was compared with the delineated MZs based on proxy variables using the Jaccard similarity index. The Jaccard index quantifies the degree of similarity between two finite sets (i.e., two MZs) and is a widely used method to compare the results of image segmentation analysis. The Jaccard index is formally defined as the size of the intersection divided by the size of the union of two finite sets: where J(A, B) is the Jaccard index between sets A and B, ∩ represents the intersection, and ∪ represents the union between the two sets. The Jaccard index ranges from 0 to 1, where a similarity of 1 means perfect matching between the two sets and a similarity of 0 represents complete discrepancy between the two sets. The analysis was conducted using the jaccard function from the Matlab Image Processing Toolbox version 10.3.

Determination of the location of SMS
To objectively define the tentative sensor location within each delineated MZ, we defined a sensor location index (SLI) that accounts for both the FCM membership grade and the distance-to-edge of the MZ for each grid cell: where SLI is the sensor location index at grid cell i, is the Euclidean distance from grid cell i to the nearest boundary of the delineated MZ, max is the maximum Euclidean distance of any grid cell in the delineated zone to the nearest boundary of the zone (i.e., distance of the centermost grid cell of the zone), and is the membership grade at the grid cell i. In other words, the SLI finds the center-most grid cell weighted by the membership grade. Both the relative distance ( ∕ max ) and membership ( ) components of the SLI are continuous and range between 0 and 1, meaning that the resulting SLI is also bounded in the range [0, 1] and is comparable across different MZs. Then, the tentative geographic coordinates for the optimal sensor location in each delineated zone can be found by searching the grid cell with the maximum SLI value (Equation 4).
The proposed method assumes that the location of single soil moisture-informed SMS can represent the soil moisture conditions of the entire MZ. The assumption rests on the notion that the zoning process substantially reduces the soil moisture spatial variability within each delineated MZ (Barker et al., 2017;Reyes et al., 2019). A distance transform was used to compute the Euclidean distance between each grid cell and the zone boundaries used the Matlab function bwdist from the Image Processing Toolbox version 10.3 (Maurer et al., 2003). We considered the distance to the edge of the MZ because of two reasons: (a) edges of agricultural fields are usually exposed to frequent traffic resulting in higher soil compaction that can affect soil water redistribution and storage, and (b) boundaries between zones constitute transitional zones that may not fully represent either zone (e.g., membership grades around a value of .5). Elevation data for each field was obtained from USGS elevation models (USGS, 2017) at 10-m spatial resolution and downscaled to a 5-m grid using a cubic spline interpolation method to match the field grid. In previous stages of this project, we considered using a digital elevation map to compute the topographic wetness index, which is typically used for representing wetness conditions at the watershed level. However, the relatively small elevation gradients and limited topographical extent to estimate the upslope area within the boundaries of the agricultural fields resulted in ill-conditioned maps that eventually caused us to favor the use of elevation maps over the topographic wetness index.

Determination of soil physical properties
In order to characterize the fraction of sand, silt, and clay across the agricultural fields, we collected a total of 113, 91, and 106 disturbed soil samples from the top 12 cm across Fields A, B, and C, respectively (Figure 1). Disturbed soil samples were collected using a shovel and then homogenized in a plastic container before bagging the samples. Soil samples were oven dried at 105˚C for 48 h and ground to pass a 2-mm sieve. Particle size analysis was conducted using the hydrometer method (Gavlak et al., 2003). To characterize the soil water retention properties of each MZ, we collected three undisturbed soil samples in the 0-to-5-cm layer per delineated zone (total of 18 samples across all three fields) using 250-cm 3 stainless-steel rings. The laboratory procedure for measuring soil water retention curves consisted of first saturating each soil sample in a 5 mM CaCl 2 solution.

Vadose Zone Journal
F I G U R E 1 Soil textural classes of the top 12 cm. Black crosses (+) represent soil sampling location across each field. In total, our study spanned five major soil textural classes (excluding a single observation of sandy clay loam soils) After reaching saturation, each sample was mounted on a sensor head equipped with two precision pressure transducers and two mini-tensiometers (2.5 cm and 5 cm long; Hyprop 2, Meter Group) connected to a desktop computer. The soil matric potential of each mini-tensiometer was recorded using the Hyprop-View software (Meter Group) while the samples remained exposed to evaporation at laboratory ambient conditions (temperature = 23 ± 2˚C, relative humidity = 21 ± 5%) (Schindler et al., 2010;Wind, 1966). Sample mass was recorded three times per day during the first 3 d of the sample dry-down and twice a day during the remaining period. Measurements with the Hyprop 2 were terminated after the longer mini-tensiometer (i.e., closest to the sample surface) reached its cavitation point (i.e., ∼80 kPa). The next step consisted of dismounting the samples from the Hyprop 2 and then extracting small soil subsamples along the wetting front of the sample to determine soil water potential using a dewpoint soil water potential meter (WP4C, Meter Group). Soil samples were weighed individually and oven dried to obtain the gravimetric water content of each sample. Data were exported and a soil water retention model (van Genuchten, 1980) was fitted using the lsqcurvefit function in Matlab 2018b (Mathworkss). To quantify the osmotic potential of the soil solution, six samples per field were measured using the soil water potential meter. The soil solution was obtained from centrifugation of 5 g of ground and oven-dried soil and 5 ml of deionized water into a 10-ml centrifuge tube for 30 min at 2,500 rpm. Similar to previous studies (Bittelli & Flury, 2009), the osmotic potential of all samples was below the detection limit of 0.01 MPa of the calibrated dewpoint soil water potential meter, so we assumed that the readings of soil water potential were approximately equal to the soil matric potential.

Field-level soil moisture spatial variability
Field-average soil moisture across all surveys ranged from 0.247 to 0.388 m 3 m −3 , with Field B exhibiting drier average soil moisture conditions than Fields A and C, with mean values ranging between 0.247 and 0.314 m 3 m −3 (Table 2). These drier surveys in Field B were likely associated to the predominance of sandy loam and loam soils (Table 3) with low soil water-holding capacity. On the other hand, Field C had wetter field-average soil moisture conditions than Fields A and B, with mean survey values ranging from 0.323 to 0.388 m 3 m −3 (Table 2). Wetter conditions in Field C were likely associated with higher clay (26%) and silt (51%) contents (Table 3). Across all surveys, the CV ranged between 6.8 and 20.3%, showing a relatively strong (r 2 = .76, Table 2) inverse linear relationship with field-average soil moisture. Our results are consistent with previous studies across multiple geographic regions and land covers (Brocca et al., 2007;Famiglietti et al., 1998;Hedley & Yule, 2009;Tague et al., 2010;Western et al., 2003), in which the soil moisture spatial variability was greatest at intermediate soil moisture conditions (∼0.25-0.30 m 3 m −3 ) and lower at near saturation conditions. In our study, the highest CV of 20.3% corresponded to the field survey on 5 June 2019 in Field B (Table 2, Figure 2b), with a field-average soil moisture of 0.247 m 3 m −3 . Despite this survey being the driest in our dataset, a value of 0.247 m 3 m −3 typically represents intermediate soil moisture conditions in these type of fine-textured soils. One reason for the lack of surveys with low average soil moisture conditions at the field level is attributed to the high antecedent precipitation, with 9 out of 16 surveys F I G U R E 2 Histograms (a-c) and fitted semivariogram models of the driest (d-g) and wettest (h-j) soil moisture surveys in each field. Mean (±1 SD) values are provided for each survey, and the mean soil moisture of each survey is also indicated with a vertical dashed line. Orange histogram represents the driest survey and blue histogram represents the wettest survey having more than 90 mm of rainfall in the preceding 30 d of the survey and only three surveys with zero precipitation in the preceding week of the survey (Table 2).
At intermediate soil moisture conditions, the soil moisture spatial variability is known to be jointly dominated by microtopography and particle size distribution (Famiglietti et al., 1998). On the other hand, the lowest CV of 6.8% was observed in the wettest field survey on 8 Apr. 2019 in Field C (Table 2, Figure 2c), where soil moisture spatial variability is known to be dominated by soil physical properties such as soil porosity and hydraulic conductivity (Famiglietti et al., 1998). Grayson et al. (1997) argued that soil moisture patterns are expected to be highly organized under wet conditions, whereas in drier conditions the soil moisture patterns tend to be more irregular and influenced by actively growing vegetation, soil properties, and microtopography.
The soil moisture spatial dependence was analyzed by fitting a semivariogram to all points in each survey. Soil moisture surveys exhibited a nugget-to-sill ratio < 25%, which can be considered as a strongly spatially dependent process dominated by intrinsic properties such as soil texture and mineralogy (Cambardella et al., 1994). The median level of spa-tial autocorrelation was concentrated at distances of 180 m for Field A, 150 m for Field B, and 280 m for Field C (for reference, all fields are approximately 800 m long in the north-south direction). The shortest range value of 77 m was observed on 8 Apr. 2019, and the largest range value of 325 m was observed on 1 Apr. 2019, both in Field C and in surveys only seven days apart (Table 2). Between these two surveys, a 13-mm precipitation event shifted the field-average surface soil moisture from 0.347 to 0.388 m 3 m −3 , which reduced the CV by ∼40% and the semivariogram range parameter by a factor of four compared with those metrics on 1 Apr. 2019 (Table 2). It is worth noting that because soil moisture observations in different surveys were made in different positions of the field, the spatial variation of the change in water content is superimposed on the variation in the water content on each measurement occasion. Experimental semivariograms of the driest surveys (Figure 2d, f, g) and the wettest surveys (Figure 2h, i, j) showed a relatively small nugget value, which is expected for dense sampling campaigns with a separation distance between observations approaching a value of zero meters. Despite the wide gamut of semivariogram range values, we did not find any meaningful correlation between semivariogram range and mean volumetric water content. Our findings suggest that the range of soil moisture spatial autocorrelation for agricultural fields at intermediate soil moisture conditions may be dictated by other factors such as soil physical properties and topography rather than soil moisture perse.

Delineation of MZs using soil moisture
The use of a recursive silhouette method coupled with the FCM consistently divided the three studied fields into two soil moisture-based MZs that captured the most salient soil moisture regimes across fields with different spatial extents, soil texture, and topographic conditions (Figure 3). Previous studies using the silhouette method in combination with clustering analysis also found that the optimal number of field MZs ranged between two to three based on multiple soil, crop, and topographic variables (Reyes et al., 2019). Maximum silhouette criterion scores were between 0.60 and 0.68 (Figure 4), indicating well-structured clusters (Kaufman & Rousseeuw, 1990) that effectively captured the distinct soil moisture regimes of each field. Fields A and B were divided into two MZs of even area, but Field C was divided into a larger wet MZ spanning nearly 76% of the field area (Table 3). This uneven distribution in Field C is attributed to the localized presence of loam soils in the lowland area of the field that resulted in drier soil moisture conditions, which, as mentioned earlier, likely contributed to the high kurtosis in the soil moisture distributions of Field C. As expected, the clustering analysis resulted in MZs with clearly distinct average soil moisture conditions and reduced intrazone soil moisture spatial variability compared with that of the entire field (Whelan & McBratney, 2000). The dry MZs presented average soil moisture conditions ∼8% lower, whereas the wet MZs showed between 4.4 (Field C) and 9% (Field B) higher soil moisture compared with the average of the field. The intrazone spatial variability quantified by the SD was substantially reduced across all three fields by about 50%. For instance, in Field C, the soil moisture SD of the dry MZ was 0.026 m 3 m −3 and the SD of the wet zone was 0.017 m 3 m −3 , values that are 55% of the soil moisture SD of 0.039 m 3 m −3 at the field level. Despite generating new more homogeneous MZs, a distinct feature of the FCM is the ability to retain information about intrazone spatial variability based on the assigned membership grade. In contrast, hard clustering techniques (i.e., K-means) lose important information about the intrazone spatial variability and cannot reveal less representative or transitioning areas between MZs. Although several studies have reported satisfactory results using K-means clustering (Arnó et al., 2012;Cambouris et al., 2006;Haghverdi et al., 2015), fuzzy algorithms are often regarded as a more F I G U R E 3 Delineated management zones (MZs) using the fuzzy c-means (FCM) clustering technique based on surface (0-12 cm) soil moisture observations for each field. The membership grade (1 is highest and 0 is lowest membership) indicates the degree of belonging of each grid cell to a specific MZ. Triangle markers represent the location of undisturbed soil samples for characterizing the soil water retention properties of each MZ F I G U R E 4 Silhouette values used to define the optimum number of management zones at each field. Silhouette values range between −1 and 1, where the highest value indicates the optimal number of clusters appropriate approach to classify for continuously variable natural phenomena (Burrough, 1989).
To further investigate the controlling factors of each MZ, we also collected undisturbed soil samples to characterize the soil water retention properties ( Figure 5). In general, wet MZs were characterized by fine-texture soils (e.g., clay loam, silty clay loam) presenting a 14% higher effective water content at saturation, 38% higher soil water content at field capacity (−10 kPa), and 207% higher soil moisture content at the permanent wilting point (−1,500 kPa) compared with dry MZs characterized by coarser soil textures (e.g., sandy loam and loam soils; Table 3). It is worth noting that in Field B, the wetter MZ was characterized by loam soils and the dry MZ was characterized by even drier soil moisture conditions due to the presence of a coarse sandy loam soil.

Delineation of MZs using proxy variables
Because intensive soil moisture surveys are time consuming and labor intensive, we also analyzed MZs delineated only using the two surveys with most contrasting field-average soil moisture content and common proxy variables such as elevation, clay fraction, and sand fraction. Regardless of the variable, this new classification using surrogate variables resulted in two contrasting MZ for the three fields. The Jaccard (J) index was calculated to analyze the spatial similarity between soil moisture-based MZ and proxy variables (Table 4). Delineated MZs based on only two soil moisture surveys presented the highest J similarity with values ranging from 0.60 (Field C, dry MZ) to 0.93 (Field B, wet MZ), and an average across all fields of 0.79. Thus, using two soil moisture surveys collected under contrasting field conditions arises as a plausible alternative to characterize the field soil moisture spatial variability and delineate MZ using a limited number of observations. For instance, intensive soil moisture surveys could be conducted before and after a major rainfall event to ensure contrasting field conditions between surveys. In Vadose Zone Journal F I G U R E 5 Soil water retention curves for each management zone. Hyprop mini tensiometers cavitated at tensions near the last observation at a matric potential (ψ m ) around −80 kPa T A B L E 4 Jaccard index values used to compare the spatial similarity between soil moisture-based management zones (MZs) with delineated MZ using proxy variables: only two soil moisture surveys, field elevation, and soil clay and sand fraction contrast, delineated MZ using elevation showed the lowest similarity values with an average of J = 0.48. In fields with small changes in topographic conditions (<1% slope; e.g., Field A), elevation seems to be weakly (J < 0.6, Table 4) associated with soil moisture spatial patterns (Figure 6d). However, in fields with well-defined topographic conditions between MZs (e.g., Field C), elevation could be used as a proxy variable to characterize soil moisture-based MZ reaching J similarity values of 0.82 (Figure 6f). Another tested alternative was surface soil clay and sand fractions. In this case, we did not use the soil textural class from public databases (e.g., Soil Survey Geodatabase, SSURGO), but rather the soil texture information derived from particle size analysis using ∼100 soil samples from each field. The resulting MZs revealed that interpolated maps of clay and sand content were able to successfully capture the main spatial pattern of soil moisture MZs with four out of six MZs exhibiting a J similarity value >0.5 (Table 4, Figure 6). Our results are in agreement with previous studies in which soil texture was a key variable to delineate MZs (Cambouris et al., 2006;Moral et al., 2010;Reyes et al., 2019).
During the fallow period in 2019, a survey of apparent electrical conductivity (EC a ) was available to us from the owner of Field C. The EC a survey was conducted with an electromagnetic inductance device (Model EM38, Geonics) that collected information for both shallow (0.5-m depth) and deep (0.7-m depth) soil layers. However, the low correlation between EC a and clay content (r 2 = .14%) and between EC a and surface volumetric water content (r 2 = .45) dissuaded us from using this single EC a survey to delineate MZs. The weak correlations in our dataset agree well with previous studies also involving observations of soil moisture on the same day of the EC a survey, in which the correlation between EC a and gravimetric water content for a single survey were also weak (r 2 < .21; García et al., 2012). Single EC a surveys have been used to characterize the spatial variability of agricultural fields (Cambouris et al., 2006;Haghverdi et al., 2015; F I G U R E 6 Alternative delineation of management zones (MZs) using two soil moisture surveys (top), elevation (middle), and combined clay and sand fraction information (bottom) for Fields A (left), B (middle), and C (right) Johnson et al., 2003;Reyes et al., 2019), but MZs based on a single EC a survey can result in high intrazone soil moisture spatial variability and in highly fragmented zones of small area that would make precision farming operations difficult (Hedley & Yule, 2009). Recent evidence shows that multiple ECa surveys are required to improve the accuracy of soil moisture estimations based on EC a surveys (García et al., 2012;Huang et al., 2017).

SMS location
The proposed index allowed us to define the tentative location within each MZ to install a SMS (Figure 7). The tentative location is represented by the geographic coordinates of the center-most grid cell that best represents the MZ based on the assigned membership grade. The approach effectively ignores grid cells of high membership grade that are located at F I G U R E 7 Optimal locations for installing one soil moisture sensor per delineated management zone based on the sensor location index (SLI). The SLI ranges from 0 to 1, where high values represent centralized grid cells within each management zone with high membership grade or near the edges of the MZ and centered grid cells of low and moderate membership grade. While deploying SMS near the edge of the field can be convenient, agricultural fields usually exhibit a notorious edge-effect that could affect the representativeness of SMS near the edges of the MZ (Carlesso et al., 2019). The edge effect is usually caused by a combination of factors, including increased soil compaction due to excessive traffic with farming equipment, higher plant population due to overlapping of planter passes, and a different microclimate relative to the bulk of the crop (Augustin et al., 2020;Carlesso et al., 2019). Some of these variables, like soil compaction, could be included in the analysis by adding observations of penetration resistance, which could help in better delineating headland zones. The computing of the SLI is simple, bounded by a specific range (i.e., 0 to 1), and comparable among MZs. Naturally, the optimal sensor location is largely governed by the variables used during the clustering process. In our case, the direct observation of surface soil moisture accounted for the salient soil moisture spatial patterns as a result of the interaction of local environmental conditions, topographic features, soil texture, and management practices (e.g., presence of residue cover). The incorporation of alternative factors partially related to soil moisture such us organic matter, crop F I G U R E 8 Sensitivity analysis for the tentative location of soil moisture sensors using all soil moisture surveys, only the driest and wettest soil moisture surveys, and proxy variables such as elevation and the fraction of clay and sand. Filled markers represent the sensor location in the dry management zone (MZ), and hollow markers represent the sensor location in the wet MZ. The circle markers are the tentative location for each zone using all soil moisture surveys and represent the benchmark location yield, apparent electrical conductivity, soil compaction, and vegetation indices could simplify the methodology, but these factors could also introduce spatial information unrelated to soil moisture, affecting the optimal location of sensors for soil moisture monitoring. To explore the impact of additional variables, we conducted a sensitivity analysis comparing the final sensor location using the proposed method base on soil moisture against common surrogate variables (Figure 8).
The first alternative consisted of using only the soil moisture surveys with the highest and lowest mean volumetric water content in each field. Not surprisingly, this alternative resulted in two well-defined MZs that exhibited similar spatial patterns to those of using the full dataset of soil moisture surveys (Table 4, Figure 5). The largest discrepancy was observed in the wet zone of Field C, in which the tentative sensor location using only two soil moisture surveys differed by 147 m relative to the benchmark location using all soil moisture surveys, and the lowest discrepancy was observed on the dry zone of Field A, where the two approaches differed by only 12 m, a relatively minor distance within the context of large production fields. A second alternative that we explored consists of delineating soil moisture MZ using observed clay and sand content. Because of the strong influence of particle size on the water holding capacity of soils, clay and sand fractions are often used as a proxy for soil moisture spatial variability. The results using soil texture showed an increased discrepancy compared with the benchmark location using all soil moisture surveys, with tentative locations ranging from 37 to 463 m. The lowest discrepancies were observed in Field C, with values below 50 m on both MZs, whereas Field A presented the largest discrepancies with values of 199 m for the dry MZ and 463 m for the wet MZ. As a third alternative, field elevation was used to determine the sensor deployment location, which resulted in similar discrepancies to those observed using clay and sand content. These results suggest that soil texture or elevation are able to capture the bulk shape of the soil moisture MZ in some scenarios, but other factors like microtopography, soil compaction, and stubble distribution may be affecting the soil moisture spatial patterns that can only be captured by direct observations of soil moisture.
The scope of our study was soil moisture spatial patterns as a result of the interplay between environmental conditions, soil physical properties, topography, and crop residue during fallow periods and early crop stages. Despite not accounting for vegetation effects on soil moisture spatial patterns, intensive soil moisture surveys in periods with no or little vegetation revealed areas of the field with clearly distinct soil moisture regimes. Because of the strong link and control of soil texture (Cosh & Brutsaert, 1999) and topography on the soil moisture spatial patterns, we speculate that it is unlikely that the presence of vegetation would alter the number and total area of the delineated soil moisture-based MZs, particularly in the absence of severe crop nutritional deficiencies or diseases, but more research is needed in this area to confirm our assumptions during the growing season. The presence of an actively growing crop could affect the membership value of the grid cells, the boundaries of the MZs, and consequently the sensor deployment locations. Another limitation of our dataset is the use of near-surface (i.e., 12 cm) soil moisture observations rather than rootzone soil moisture observations, an issue that was entirely related to the time and labor involved in the intensive collection of in situ soil moisture across large production fields. Future studies aiming at optimizing the number and deployment location of SMS should explore profile-level soil moisture dynamics and the impact of vegetation water uptake on the soil moisture spatial patterns at deeper layers. Because of the substantial amount of time and labor involved in intensive field surveys, future studies could focus on comparing in situ soil moisture observations with remote and proximal multispectral images of bare soil and vegetation, multiyear yield monitor data, time-lapse surveys of apparent electrical conductivity, and cosmic-ray neutron rover surveys of field soil moisture (Schrön et al., 2018). Physically based models (e.g., HYDRUS; Šimůnek et al., 2008) in combination with in situ observations can span a wider range of soil moisture scenarios and also arise as potential strategy for identifying the minimum number and tentative location of SMS.

CONCLUSIONS
The combined use of the FCM unsupervised clustering technique coupled with the silhouette cluster evaluation method allowed us to objectively identify and delineate a specific number of field MZs based on intensive surface soil moisture observations across three agricultural fields of varied soil texture and topographic conditions. Across the three studied sites, two soil moisture MZs were sufficient to capture the dominant drier and wetter soil moisture regimes mainly dictated by soil texture, landscape position, and soil hydraulic properties. A sensor location index that combined the FCM membership grade and the distance-to-edge of the MZ resulted effective to define the tentative location of SMS within the previously delineated soil moisture-based MZs. Intensive surveys of surface soil moisture observations using a portable soil water reflectometer revealed the complex spatial patterns emerging from the interplay of weather, topography, management, and soil physical properties. Although MZs delineated using the percentage of clay and sand showed better agreement compared with using elevation alone, neither the use of soil texture information nor field elevation were able to closely match the tentative location for the deployment of SMS compared with the benchmark using multiple in situ soil moisture surveys. Our study provides evidence that proxy variables may not accurately represent the underlying field soil moisture patterns and that, when possible, in situ soil moisture observations should be favored, or at least included, to guide the delineation of field MZs for improved zoning of the water management and for guiding the optimal location of a limited number of SMS.

C O N F L I C T O F I N T E R E S T
The authors declare no conflict of interest.