On the Radiative Transfer Model for Soil Moisture across Space, Time and Hydro-Climates

A framework is proposed for understanding the efficacy of the microwave radiative transfer model (RTM) of soil moisture with different support scales, seasonality (time), hydroclimates, and aggregation (scaling) methods. In this paper, the sensitivity of brightness temperature TB (H- and V-polarization) to physical variables (soil moisture, soil texture, surface roughness, surface temperature, and vegetation characteristics) is studied. Our results indicate that the sensitivity of brightness temperature (V- or H-polarization) is determined by the upscaling method and heterogeneity observed in the physical variables. Under higher heterogeneity, the TB sensitivity to vegetation and roughness followed a logarithmic function with an increasing support scale, while an exponential function is followed under lower heterogeneity. Surface temperature always followed an exponential function under all conditions. The sensitivity of TB at H- or V- polarization to soil and vegetation characteristics varied with the spatial scale (extent and support) and the amount of biomass observed. Thus, choosing an H- or V-polarization algorithm for soil moisture retrieval is a tradeoff between support scales, and land surface heterogeneity. For largely undisturbed natural environments such as SGP’97 and SMEX04, the sensitivity of TB to variables remains nearly uniform and is not influenced by extent, support scales, or an upscaling method. On the contrary, for anthropogenically-manipulated environments such as SMEX02 and SMAPVEX12, the sensitivity to variables is highly influenced by the distribution of land surface heterogeneity and upscaling methods.


Introduction
The remote sensing community aims to develop soil moisture products at various resolutions, as soil moisture finds application from local (e.g., crop water management) and regional (e.g., flood and drought forecasting) to global scales (e.g., meteorology, climate dynamics) [1]. The active or passive sensors used for estimating soil moisture come with their own challenges such as accuracy, Radio Frequency Interference (RFI), high sensitivity to vegetation/roughness, spatio-temporal coverage, big data challenges, etc., [2,3]. Developing fine resolution soil moisture products can be achieved either by downscaling retrieved coarse-resolution soil moisture or downscaling the observed brightness temperature (T B ) followed by soil moisture retrieval [4][5][6][7]. Several scaling algorithms have been developed in the past using data fusion techniques integrating information from various scales (point/airborne/satellite), platforms (MODIS/LANDSAT/AVHRR etc.), sensors (active/passive), frequencies (P, L, C, and X), and land surface variables (surface temperature, NDVI, etc.) [8][9][10][11][12]. For most cases, coarse-scale knowledge is estimated from smaller spatial scale features in which upscaling has reduced to the problem of a change of support scale. Western and Bloschl (1999) [13] defined the scale

Heterogeneity Observed in Various Hydroclimates
The effects of varying heterogeneity on brightness temperature are studied in four different hydroclimates using data from field campaigns, namely SGP97 (Southern Great Plains'1997), SMEX02 (Soil Moisture Experiment'2002), SMEX04 (Soil Moisture Experiment'2004), and SMAPVEX12 (Soil Moisture Active Passive Experiment'2012). These field campaigns were conducted over approximately a 4-8 weeks' window, primarily to validate soil moisture retrieval algorithms over a wide range of soil and vegetation conditions [31][32][33]. The airborne remote sensing data for these field campaigns were observed at either 800 m or 1.5 km resolution, which is the base resolution for our analysis. The 800 m pixels are upscaled to support scales of 1.6 km, 3.2 km, 6.4 km, and 12.8 km, and 1.5 km pixels are upscaled to 3 km and 9 km to coincide with two of the Soil Moisture Active Passive (SMAP) target products. The variables such as soil moisture, surface temperature, vegetation water content, and soil texture are obtained from various sources, which are discussed under each field campaign. Variables such as surface roughness, vegetation structure, and single scattering albedo, which are available only at sparse locations, are assumed to be scalars with a uniform probability distribution whose ranges are estimated from extensive field measurements.
2.1.1. Southern Great Plains Experiment'1997 (SGP'97), Oklahoma The SGP'97 hydrology experiment was conducted in the sub-humid (Köppen climate classification-Cfa) [34] region of central Oklahoma from 18 June through 17 July 1997. The topography of the region is moderately rolling and soils vary with a wide range of textures, including large regions of both coarse and fine textures, and with an average annual rainfall and temperature of 750 mm and 288.85 K, respectively. The land use is dominated by rangeland and pasture with significant areas of winter wheat and other crops. The L-band Electronically Scanned Thinned Array Radiometer (ESTAR) was used to measure soil moisture at 800 m × 800 m resolution. An 800 m resolution daily effective soil temperature maps were generated using the ground-based soil temperatures and grid-based resampling program. The vegetation water content (VWC) map at 800 m resolution was developed from normalized difference vegetation index (NDVI) values computed using the Landsat TM imagery collected on 25 July 1997 [35]. A combination of soil geographic database (STATSGO) and soil texture samples were used to develop a resampled texture map at 800 m resolution [23].

Soil Moisture Experiment'2002 (SMEX02), Iowa
The SMEX02 experiment was conducted in Walnut Creek watershed, Iowa, from 23 June through 12 July 2002. This region is classified as a humid climate (Köppen climate classification-Dfa) with an average annual rainfall of 835 mm and a temperature of 283 K. A considerable amount of variability in soil texture is observed in the region, ranging from fine sandy loam to clay with the majority classified as silt loam with relatively low permeability [36]. Corn and soybean crops are the two dominant land covers of Iowa, covering approximately 50% and 40%, respectively. The C-band Polarimeteric Scanning Radiometer (PSR-C) retrieved soil moisture at 800 m × 800 m resolutions [37] and was used for the analysis. A total of four days (DOY 178, 182, 186, and 188) were selected to capture the temporal variability in land surface heterogeneity due to crop growth, temperature, and soil moisture changes. The Soil Survey Geographic Database (SSURGO) and Landsat-7/5 E/TM derived vegetation water content at 30 m that was resampled to 800 m × 800 m grids. The averages of the surface and 1-cm depth soil temperature were interpolated through kriging to obtain coverage for a grid size of 800 m.

Soil Moisture Experiment'2004 (SMEX04), Arizona
The SMEX04 field campaign was conducted with a primary focus on Walnut Gulch (WC) Experimental Watershed near Tombstone, which is operated by the Agriculture Research Service (ARS), U.S. Department of Agriculture (USDA), between 2 August and 27 August 2004 [31]. The climate of this region is classified as semi-arid (Köppen climate classification-BSh), with a mean annual temperature of 290.85 K and receiving an average annual precipitation of 350 mm. The soils are generally well drained calcareous and gravelly loam with large percentages of rock ranging from nearly 0% on shallow slopes to over 70% on the very steep slopes. Vegetation of this region is mainly comprised of two-thirds of shrubs and the remaining one-third is grassland. The polarimetric scanning radiometer (PSR) C-band retrieved soil moisture and regularly spaced grid of 800 m resolution was used. The ground-based temperatures (surface and 5 cm temperatures are averaged) and MODIS (MOD11A1) surface temperature product are interpolated (kriging) to generate an 800 m resolution surface temperature map. The 800 m resolution vegetation water content map was developed by interpolating high resolution Landsat TM and regularly ground sampled VWC [38]. The SSURGO data at 30 m were resampled to obtain a soil texture map at 800 m resolution. A total of four days (DOY 221, 222, 225, and 226) were selected based on the availability of datasets to study the temporal variability in land surface heterogeneity. The surface roughness data were collected extensively over the Walnut Gulch (WG) watershed and sparely over the entire SMEX04 region using a pin board. SMAPVEX12 was conducted from 6 June to 17 July 2012, in Winnipeg, Manitoba (Canada), which is classified as having extreme humid climate (Köppen climate classification-Dfb), with average an annual precipitation and temperature of 521 mm and 289 K, respectively. The region captures a wide variety of variability in land cover and soil texture within a few kilometers while moving from southeast to northwest [33]. The agricultural crops (mainly cereals, soybeans, canola, and corn) dominate the southeast with heavy clay soils, whereas the north and east see more lighter sandy and loamy soils with mixed land cover including forest and perennial pastures along with a few crop fields. The analysis was conducted for all days (DOY 164,167,169,174,177,179,181,185,187,192,195,196,and 199) using the Passive Active L-band System (PALS) derived soil moisture and surface temperature at 1.5 km resolution [39]. The VWC available at 5 m resolution estimated from the Normalized Difference Water Index (NDWI) using SPOT and RapidEye satellite overpasses was resampled to 1.5 km grids.

Soil Moisture Retrieval Algorithm
The approach followed in this study was based on forward modeling using a first-order microwave radiative transfer model (RTM), also called a tau-omega (τ-ω) model [40], which is described in Equation (1) to simulate brightness temperature (T B ), where T B(p, f ,θ) is the brightness temperature [K]; T e f f is the effective surface temperature [K]; T c is the effective vegetation temperature [K]; e p,θ, f is the emissivity of the (rough) soil surface; r p, f ,θ is the rough surface reflectivity; τ p, f is the nadir optical depth; ω p, f ,θ is the single scattering albedo. The subscripts p, θ, and f denote the polarization, angle of incidence, and frequency of measurement, respectively. It is assumed that the canopy and soil temperature are equal during morning thermal crossover, which also corresponds to the planned time of observation for SMAP at around 6:00 h local sun time.
The microwave signal incident emitted from the soil-vegetation layer is modified, scattered, and attenuated due to the physical and structural properties of this medium. The radiative transfer (Equation (1)) is a first-order approximation and summation of three components: (i) the direct emission by soil and one-way attenuation by canopy (the first term), (ii) direct upward emission by canopies (the second term), and (iii) emission by plants and reflected by soil and thereafter attenuated by vegetation (the third term). The physical attributes of soil layer that influence the microwave emissions are surface reflectivity r p, f ,θ and effective soil temperature T e f f , as given in Equation (1). The physical attributes of the vegetation that influence the microwave emissions are: the optical thickness of the vegetation τ p, f , and the single scattering albedo ω p, f ,θ . The formulation in Equation (1) is based on the assumptions that, (i) the diffuse scattering is ignored; (ii) the air-vegetation reflectivity is assumed zero, thus ignoring losses at the boundary; and (iii) the refractive index of the vegetation is approximately equal to air, which allows us to use the air-soil reflectivity in Equation (1). The extinction optical thickness of the vegetation τ p is estimated as a product of vegetation water content (VWC) and the structure of the canopy (b). Due to the complex geometry of the natural canopies, approximate values are used for the canopy structure (b). The vegetation water content is estimated using the Normalized Difference Vegetation Index (NDVI) and regression equation [41]. The τ p and ω p are polarization dependent because there are differences in propagation of H-polarized and V-polarized waves under different canopy structures [42][43][44][45]. The surface roughness represents the small-scale surface undulations. These surface undulations are typically measured using a pin-board or a grid board, where the variation in elevations of metal pins are considered due to the rough surface underneath. However, due to the infeasibility of doing these measurements operationally at airborne footprints (1.5 km, 3 km, etc.). Neelam et al., 2020 [46], developed a multi-scale surface Remote Sens. 2020, 12, 2645 5 of 22 roughness model based on variables such as soil moisture, leaf area index (LAI), clay fraction (CF), and topographic curvature. Nonetheless, the current SMAP operational algorithm uses static surface roughness parameters based on IGBP land cover classes [47].

Global Spatial Sensitivity Analysis: Sobol Method
The Global Spatial Sensitivity Analysis (GSSA) [27,30] relies on Sobol methods [48], which can deal with non-linear and non-monotonic relationships between inputs and output [49][50][51][52]. These methods are based on the functional decomposition of variance (ANOVA) of the model prediction (Y) into partial variances caused due to each model input (X) (either considered singularly or in combination), i.e., where V i (Y) is the share of the output variance explained by the ith model input, V ij (Y) is the share of the output variance explained by the interaction of the ith and jth inputs, and k is the total number of inputs.
The partial variances are explained as , and so on for higher order interactions, where E and V are the expectation and variance operators; refer to the Appendix A for more details. The so-called "Sobol" indices or "variance-based sensitivity indices" are obtained as follows: , the amount of variance of Y explained by the interaction of the factors X i and X j (i.e., sensitivity to X i and X j not expressed in V i nor V j ); Total Sensitivity Index S Ti = S i + i< j S ij + i<j<l S ijl + . . . S 1,2,...k , which accounts for all the contributions to the output variation due to factor X i (i.e., the first-order index plus all its interactions). In other words, the S Ti index is defined as a summation of the main, second, and higher order effects, which involves the evaluation over a full range of parameter space. If S i and S Ti are equal and the sum of S i (and thus S Ti ) is 1, then the model is additive (linear) in nature, otherwise if S Ti is greater than S i and S i < 1 or S Ti > 1, then the model exhibits non-linearity. Thus, linearity or non-linearity of model through time, scale, and hydroclimate can be determined from a summation of S Ti s indices.
Let Y = f (X) be a spatial model, then Y is the model output and X are random independent model inputs, where X is composed of U and W, with U defined as random vectors and W defined as two-dimensional spatial maps. Using the similar configuration, Y is the simulated brightness temperature, f is the radiative transfer model, and X K is the input parameter vector with K = [X 1 , X 2 , X 3 , X 4 ] as [Surface Roughness-RMS height (S), Surface Roughness-Correlation length (L), Vegetation Structure (B), Scattering Albedo (ω)] and W as spatial maps [Soil Moisture (SM), Clay Fraction (CF), Surface Temperature (TSURF), Vegetation Water Content (VWC)]. The steps used to carry out the analysis are described in Figure 1. The variables such as RMS height (S) and correlation length (L) are sampled from a uniform distribution whose ranges are estimated from field measurements. The variables such as vegetation structure (B) and scattering albedo (ω) are also sampled from a uniform distribution whose ranges are defined based on vegetation type. However, experimental data for these variables are limited and are obtained from past literature [42,53]. The range of variability in surface roughness, vegetation parameter, and scattering albedo are assumed to remain constant at all scales, as there are no means to ascertain how those vary with scale. The mean values of land surface variables used in the analysis for various field campaigns are represented in Table 1. The spatial input variable at each support scale is defined as a real random variable with a uniform range for the number of pixels developed at that scale. That is, the numbers of pixels 'n' at each support scale are considered as equiprobable, with each pixel labeled with a unique integer in the set {1, . . . , n}. The number of pixels decreases with increasing the support scale, for example, the numbers of pixels at a 800 m support scale are more than 4000, and at a 12.6 km support scale, there are about 50 pixels. The first order and total Sobol sensitivity indices are estimated with confidence intervals using the bootstrap Remote Sens. 2020, 12, 2645 6 of 22 technique with resampling [54]. The number of samples N = 55,000 are used for model evaluations.
Here, we discuss and present only total sensitivity indices. The current SMAP Level 2 and 3 radiometer-based soil moisture algorithm [47] uses vegetation and surface roughness parameters from the look up table by IGBP class, assuming sub footprint-scale homogeneity for vegetation and roughness parameters. However, surface roughness parameters change with soil type, vegetation conditions, wind speed, amount of precipitation, etc. Similarly, vegetation parameters vary with growing season, vegetation type, and structure. Due to the heterogeneous and dynamic nature of these land surface variables (varying with space and time), assuming them to be fixed parameters will introduce high errors into soil moisture retrievals. Since the range for vegetation and roughness variables are obtained from the respective field campaign conducted over a few weeks, larger heterogeneities of land surface variables are possible to occur in real time, as we have considered in our analysis.

Upscaling Methods: Linear Upscaling vs. Inverse Distance Weighted (IDW) Upscaling
Scaling is an important aspect, and is studied in various disciplines including hydrology [55,56] meteorology [57], ecology [58], and geography [59]. The upscaling and validation of satellite soil moisture products are discussed in Crow [63][64][65]. The upscaling technique for any remotely sensed land surface variable depends on the heterogeneity exhibited by the variable and antenna gain function, which plays a significant role in capturing the heterogeneity. In this paper, we investigate the difference in sensitivity of the retrieval model to the land surface variables when the variables are upscaled linearly and through Inverse Distance Weighted (IDW) to various support scales (0.8 km, 1.6 km, 3.2 km, 6.4 km, and 12.8 km). This is important because of three main reasons: (1) as the support scale changes (upscaling or downscaling), the heterogeneity of the landscape also varies. Thus, the processes that appear linear may become non-linear and vice versa. The variables that are important at one scale and hydroclimate may become trivial at another scale and hydroclimate. (2) The non-linearity of the radiative transfer model varies with heterogeneity observed for the retrieval scale. The non-linear retrieval model may show no scale effect when the scale of analysis (support) is homogeneous [66,67], and (3) the upscaling techniques become significant in heterogeneous landscapes to represent the non-linear effects [36,[68][69][70][71]. In the case of a homogeneous environment, the scaling methods do not significantly influence the output, as is also validated in Wu and Liang Li (2009) [72]. The use of an unsuitable upscaling method or inappropriate representation of spatial distribution can lead to potentially wrong decisions. An inappropriate upscaling can have an even more profound impact if the result is used as an input for simulations, where a small error or distortion can cause models to produce false estimation.
IDW is a deterministic non-linear interpolation technique that uses the weighted upscaling of the attribute values from the nearby sample sub-grid to estimate the magnitude of the attribute at a higher support scale. The superior performance of IDW in comparison to other scaling methods is also discussed in Weber and Englund (1994) [73], Moyeed and Papritz (2002) [74], and Babak and Deutsch (2009) [75]. The weight a particular sub-grid that is assigned in the averaging scheme depends on the distance of the sub-grid from the center of the coarser grid. It is an averaging method based on the principle that a sample closer to the prediction location has more influence on the prediction than a sample farther apart. For gridded data, the IDW approach can be realized as a-3 dB function where higher weights are given to grid cells closer to the center of the pixel than grid cells further apart. This incorporates neighboring heterogeneity, which changes with the support scale. For this analysis, the weight of the point changes with the inverse square of the distance known as Inverse Distance Squared (IDS) upscaling. This technique is also used in the SMAP L1C [76] radiometer data product to average non-uniformly all the brightness temperature data samples that fall within the grid cell.

Results and Discussion
It is observed from Figures 3-6 that the brightness temperature at V-polarization is more sensitive to soil moisture, soil texture (clay fraction), and surface temperature than H-polarization under all conditions (space and time). On the other hand, H-polarization is always observed to be more sensitive to surface roughness variables (S and L) than V-polarization under all conditions. The SMEX02 and SMEX04 field campaigns were conducted using a Polarimetric Scanning Radiometer C-band (PSR/C), whereas SGP'97 and SMAPVEX12 were conducted using a Passive Active L-band System (PALS) radiometer, as mentioned earlier. Though the temporal analysis has been conducted for each hydro-climate, for the sake of brevity, we presented results for one DOY and discuss the general temporal trends observed. Neelam and Mohanty, 2015, present a detailed analysis and discussion for temporal variability in RTM sensitivity analysis. In the following sections, we discuss the variability observed in the sensitivity of the L-band brightness temperature at V-and H-polarization under: (1) Plant Structure, (2) Spatio-Temporal Scales, and (3) Upscaling and Environmental Heterogeneity.

Plant Structure
It is well understood that vegetation attenuates the microwave emission from the soil, and adds its own radiative flux to the total emission. The extent of attenuation depends on the amount of vegetation and its structure. A crop canopy may be divided into three main components: (1) leaves, (2) stalks, and (3) heads. The attenuation of radiation may, thus, occur at two layers, the upper layer consisting of heads, and a lower layer composed of leaves and stalks. The propagation characteristics are different for different canopies depending on the emission and scattering properties of these layers. Therefore, assuming polarization-independence may not be true especially for plants where the vegetation biomass is concentrated in vertically oriented dipole-like stalks and row crops. For the zero order radiative transfer model (tau-omega model), the effects of vegetation are represented by vegetation water content (VWC), single scattering albedo (ω), and vegetation structure parameter (B), where VWC captures the amount of water content in biomass, 'B' represents the canopy architecture, and ω represents scattering within the canopy.
In canopies with vertical stalks, the attenuation will be smaller for horizontally polarized waves than vertically polarized waves except at a nadir incidence angle, where stalks are not visible and appear only as small randomly oriented disks (Figure 2). For a horizontally polarized incident wave, the electric field vector E is normal to the axis of the stalk at all incidence angles, and for a vertically polarized wave, E is parallel to stalk, thereby coupling strongly [77]. This difference will be significant mainly over vegetation that exhibits some preferential orientation such as vertical stalks noticed for tall grasses, grains, wheat, corn, maize, forest, etc. For vegetation with near horizontally oriented primary branches or broad leaf plants (soybean, bean, canola, sunflower, etc.), there is stronger attenuation for a horizontally polarized wave than vertically polarized waves. It is reasonable to assume that the absorption/scattering loss factor will be approximately polarization independent for canopy components whose sizes are much smaller than the wavelength of observation.  For example, for SMEX04 (Arizona), the vegetation is mainly shrub dominated (whitethorn acacia, creosote bush, etc.) covering two-thirds of the watershed, and grass (sideoats grama, lehmann lovegrass, etc.) occupying the rest with a total mean vegetation water content 0.087 kg/m 2 . Due to the low vegetation and primarily random vegetation structure, the sensitivity to vegetation variables is insignificant for this region temporally, at all support scales and polarizations (Figure 3). Similarly, the SGP97 region (Oklahoma), which is dominantly rangeland and has a significant occupation of winter wheat crops with a total mean vegetation water content of 0.323 kg/m 2 , shows the insignificant contribution of vegetation variables towards brightness temperature (V-and H-polarization) at all spatial scales (Figure 4). SMEX02 (Iowa) has a total mean vegetation water content of 1.9 kg/m 2 . At the 800 m scale, H-polarization was observed to be (~1-4%) more sensitive than V-polarization. With an increasing support scale, V-polarization was observed to be~2% (at 1.6 km) and~15% (at 12.8 km) more sensitive than H-polarization for vegetation water content (VWC) and vegetation structure (B) ( Figure 5). Under IDW analysis, H-polarization was always temporally more sensitive than V-polarization to vegetation water content by~2-3% and vegetation structure (B) by~3-9%, except on DOY 188, where V-and H-polarization showed the same sensitivities ( Figure 5). That is, at lower vegetation water content (DOY 178, 182) and a support scale of 800 m under linear upscaling and at all support scales under IDW where the heterogeneity in VWC is maintained, a higher sensitivity to VWC and B at H-polarization is noticed. For DOY 188, with higher vegetation water content and higher support scales under linear upscaling where the heterogeneity disappeared with support scales (moving towards extent mean value), a higher sensitivity to vegetation structure (B) at V-polarization is observed. This confirmed that there is a stronger coupling with plants exhibiting definite vertical structures, and the coupling increases with increasing vegetation water content and decreasing heterogeneity. The SMAPVEX12 region shows a wide variety of agricultural crops such as cereals, soybeans, canola, corn, etc., with a total mean vegetation content water of 1.4 kg/m 2 . Unlike in SMEX02, the sensitivity to vegetation structure (B) was always observed to be higher for H-polarization at all scales. We hypothesize that the amount of vegetation water content determines the sensitivity of H-and V-polarization to the B variable. To test the hypothesis, the analysis was conducted including pixels with higher vegetation water content than those observed within the SMAPVEX12 extent, which resulted in higher sensitivity to vegetation structure (B) at V-polarization. This may occur because with increasing VWC, the dielectric constant of vegetation also increases, thereby increasing the sensitivity to V-polarization due to stronger coupling, as mentioned earlier. Thus, vegetation variable 'B' is a complex function of plant geometry and vegetation water content.

Semi-Arid (SMEX04) Hydroclimate
For semi-arid regions with low vegetation and lower evapotranspiration (lower latent heat flux and higher sensible heat flux), an increase in air and surface temperatures are observed, resulting in higher sensitivity to surface temperature. However, the heterogeneity in temperature is lost with upscaling due to lumped fluxes resulting in exponential (R 2 = 0.99) decrease (increased) in sensitivity to surface temperature (RMS height, S) with support scale. A general decrease in R 2 is noticed for the for H-pol, and this is observed across all hydroclimates, as noted in Figure 3. Due to the rocky topography of SMEX04 region, the observed surface roughness is highly random, emphasizing the higher sensitivity to RMS height (S), which increased from~11% (~63%) to~38% (~87%) for V-polarization (H-polarization), and correlation length (L) increased from~1% (~7%) to~5% (~10%). This is in contrast to SMAPVEX12 and SMEX02, where the surface roughness is more correlated due to regular agricultural practices resulting in higher sensitivity to correlation length (L) as well. The sensitivity to vegetation variables (VWC, B, and ω) and clay fractions are negligible across both polarizations (H and V), upscaling methods, and spatio-temporal scales.

Sub-Humid (SGP97) Hydroclimate
For a more natural landscape with low density biomass, as observed for SGP97, soil moisture is the most sensitive variable for V-polarization (H-polarization)~96% (~83%) at all support scales under both linear and IDW upscaling. Due to the lower variability of roughness in SGP97 and lower biomass content, incorporating neighboring pixels via IDW or using linear upscaling did not influence the sensitivity measures with support scales. The sensitivity to land surface variables (clay, surface temperature, and scattering albedo) are insignificant (~1%) at V-polarization, whereas the sensitivity to surface roughness variables (S and L) at H-polarization are observed at~8% each, and sensitivity to vegetation variables (VWC and B) are observed at~2% each (Figure 4). To further explore the influence of varying extent, the analysis is repeated for a smaller extent over the Little Washita (LW) watershed. To maintain statistically significant data points, the upscaling is conducted only up to a 3.2 km support scale. Results indicated soil moisture (in both V-and H-polarization) to still be the most sensitive variable at all scales; however, its magnitude has decreased slightly to~94% for Vpolarization and~80% for H-polarization. Sensitivity to the vegetation structure, (B) increased to~4% and~2% for H-and V-polarization, respectively, at a smaller extent than the sensitivities observed at a larger extent. This could be due to the increased attenuation caused by higher biomass content (~0.41 kg/m 2 ) for a smaller extent than observed at a larger extent (SM~0.14 v/v and VWC~0.32 kg/m 2 ).

Humid-Dfa (SMEX02) Hydroclimate
For SMEX02, the analysis is conducted for DOY 178, 182, 186, and 188, and the soil moisture sensitivity at 800 m followed a similar field scale pattern as observed in Neelam and Mohanty (2015) [26]. Unlike SMEX04 and SGP97, the mean sensitivity to soil moisture across sampling days is~28% and 13% for V-and H-polarization, respectively, at 800 m, which decreased to~8% and~1% at 6.4 km, and increased by~1-2% at 12.8 km ( Figure 5). Among other land surface variables, sensitivity to surface temperature followed a decreasing exponential function (R 2 =~0.95) and nearly disappeared at 1.6 km for H-and at 3.2 km for V-polarization. On the other hand, the sensitivity to surface roughness variables (S and L) followed a logarithmic function (R 2 =~0.97) with an increasing support scale. Similarly, vegetation water content (VWC) and vegetation structure (B) followed a decreasing (R 2 = 0.96) and increasing logarithmic function. This is because of the greater variability in surface roughness and vegetation variables observed, particularly for the anthropogenically-modified region such as SMEX02. Thus, the impact of these variables is observed longer with a support scale following a logarithmic function, compared to the surface temperature that followed an exponential function. Under IDW upscaling, the sensitivity to land surface variables with support scales changed by~1-2%. Sensitivity to variables such as soil moisture, surface temperature, and surface roughness remained uniform until 6.4 km, which later decreased by~2% for soil moisture and surface temperature, and increased for surface roughness by~2% at 12.8 km. The uniformity in sensitivity to these variables until 6.4 km under IDW is similar to the reasons explained above. However, for a selected extent, the variability in some land surface variables attains saturation at a certain scale (6.4 km), resulting in a decrease in sensitivity to those land surface variables upon further upscaling. The support scale at which this saturation is attained is dependent on the heterogeneity and extent of the analysis. The sensitivity to vegetation water content remained uniform until 3.2 km and decreased~1% with an increasing support scale, and sensitivity to vegetation structure (B) increased linearly (R 2 =0.98) by~6% with the scale. The sensitivity to clay fraction and albedo remained negligible under IDW and linear upscaling. Except for DOY 178, which was a relatively dry day, the sensitivity to clay fraction increased to~3% at V-polarization and~1% at H-polarization until 3.2 km.
It is important to note that we discuss only the total sensitivity indices (TSI) of the land surface variables. The increased sensitivity of a variable may also be due to increased interaction effects such as (S,B), (S,VWC), (S,SM), (L,VWC), and (SM,B) with growing vegetation (increased interception and scattering) [26]. This may also be one of the contributing factors of the high sensitivity to surface roughness (S and L) and vegetation structure (B). The sensitivity to scattering albedo is observed to be negligible at all scales.

Semi-Arid (SMEX04) Hydroclimate
For semi-arid regions with low vegetation and lower evapotranspiration (lower latent heat flux and higher sensible heat flux), an increase in air and surface temperatures are observed, resulting in higher sensitivity to surface temperature. However, the heterogeneity in temperature is lost with upscaling due to lumped fluxes resulting in exponential (R 2 = 0.99) decrease (increased) in sensitivity to surface temperature (RMS height, S) with support scale. A general decrease in R 2 is noticed for the for H-pol, and this is observed across all hydroclimates, as noted in Figure 3. Due to the rocky topography of SMEX04 region, the observed surface roughness is highly random, emphasizing the higher sensitivity to RMS height (S), which increased from ~11% (~63%) to ~38% (~87%) for Vpolarization (H-polarization), and correlation length (L) increased from ~1% (~7%) to ~5% (~10%). This is in contrast to SMAPVEX12 and SMEX02, where the surface roughness is more correlated due to regular agricultural practices resulting in higher sensitivity to correlation length (L) as well. The sensitivity to vegetation variables (VWC, B, and ω) and clay fractions are negligible across both polarizations (H and V), upscaling methods, and spatio-temporal scales.

Sub-Humid (SGP97) Hydroclimate
For a more natural landscape with low density biomass, as observed for SGP97, soil moisture is

Humid-Dfb (SMAPVEX12) Hydroclimate
The SMAPVEX12 region observed lower biomass content and a wide heterogeneity in soil moisture due to a more gradual variability in sensitivity to land surface variables being noticed, unlike SMEX02. To maintain statistically significant sample points, the upscaling is conducted from 1.5 km to 9 km only. The sensitivity of V-polarization (H-polarization) to soil moisture showed~74% (~52%) at 1.5 km, which decreased to~61% (~35%) at 9 km ( Figure 6). The sensitivity to clay fraction remained uniform (~6%) spatio-temporally for both V-and H-polarization. A relatively high sensitivity to clay in SMAPVEX12, unlike other regions, is due to a wide heterogeneity (~5% to~65%) in clay fraction observed in this region. The sensitivity to surface roughness variables (S and L) is uniform at~2-4% at all scales, unlike in SMEX02, where they increased with scale. This may also be due to higher surface roughness measured for SMEX02 than SMAPVEX12. The sensitivity to vegetation structure and vegetation water content increased exponentially (R 2 = 0.99) with scale. It is noteworthy that the sensitivity to VWC slightly increased with scale.
1.5 km to 9 km only. The sensitivity of V-polarization (H-polarization) to soil moisture showed ~74% (~52%) at 1.5 km, which decreased to ~61% (~35%) at 9 km [ Figure 6]. The sensitivity to clay fraction remained uniform (~6%) spatio-temporally for both V-and H-polarization. A relatively high sensitivity to clay in SMAPVEX12, unlike other regions, is due to a wide heterogeneity (~5% to ~65%) in clay fraction observed in this region. The sensitivity to surface roughness variables (S and L) is uniform at ~2-4% at all scales, unlike in SMEX02, where they increased with scale. This may also be due to higher surface roughness measured for SMEX02 than SMAPVEX12. The sensitivity to vegetation structure and vegetation water content increased exponentially (R 2 = 0.99) with scale. It is noteworthy that the sensitivity to VWC slightly increased with scale.

Upscaling and Environmental Heterogeneity
The largest up-scaling effects are observed for scenes with high heterogeneity, i.e., with strong gradients in land surface variables. The differences between the two (linear and IDW) upscaling methods are, therefore, more prominent in heterogeneous environments. A high correlation is observed between spatial heterogeneity in soil moisture and density of vegetation [23,44]. Thus, we classify the hydroclimates as heterogeneous and homogeneous environments based on biomass

Upscaling and Environmental Heterogeneity
The largest up-scaling effects are observed for scenes with high heterogeneity, i.e., with strong gradients in land surface variables. The differences between the two (linear and IDW) upscaling methods are, therefore, more prominent in heterogeneous environments. A high correlation is observed between spatial heterogeneity in soil moisture and density of vegetation [23,44]. Thus, we classify the hydroclimates as heterogeneous and homogeneous environments based on biomass density. This also complements Mo et al. (1982) [40] and Wigneron et al. (1995) [78], who suggested that the total amount of biomass within the pixel is important in affecting the microwave response.

Homogeneous Environment (Sub-Humid and Semi-Arid)
SGP'97 (Sub-Humid) and SMEX04 (Semi-Arid) are more natural landscapes composed mainly of grasslands and shrubs, with low biomass content. For SGP97, a very small percentage of higher order interactions are observed,~1% for V-polarization and~3% for H-polarization. SMEX04 observed higher order interactions of~5% for V-polarization and~9% for H-polarization. Due to the rough topography of SMEX04, a higher sensitivity to roughness variables (S and L) for H-polarization is observed, which also led to increased higher order interactions than SGP97. Because of higher sensitivity to soil moisture and nearly linear (additive) behavior of the radiative transfer model, these regions can be classified as homogenous environments. Due to the homogeneous environment, we do not observe large variabilities in sensitivity measures between linear and non-linear (IDW) upscaling methods (Figures 3 and 4). This allows the land surface variables to be scaled linearly, which is also consistent with studies of Drusch et al. (1999) [68] and Crow et al. (2000) [79]. Under homogeneous environments, SMEX04 and SGP97 can further be classified as energy rich and water rich environments, respectively. SMEX04 that observed a total mean soil moisture of~0.07 v/v and higher surface temperatures of~319 K is classified as the energy rich environment, and SGP97 that observed a relatively higher total mean soil moisture of~0.15 v/v and lower surface temperatures of~298 K can be classified as a water rich environment. For a homogeneous water rich environment, the sensitivity to land surface variables did not vary with either the support scale or the upscaling methods, and varied only by small percentages due to polarization differences. The sensitivity to land surface variables under homogeneous energy rich environments followed an exponential function with support scale and attained uniformity bỹ 3.2 km, emphasizing the linearity/homogeneity of the region.

Heterogeneous Environment (Humid Dfa and Humid Dfb)
SMEX02 (Humid Dfa) and SMAPVEX12 (Humid Dfb) are agricultural landscapes with a wide range of biomass content and observe tillage practices on a regular basis. In these regions, the sensitivity to land surface variables are complex and exhibited considerable higher order interactions (non-linearities). SMAPVEX12 observed higher order interactions of~6% for V-polarization and~10% for H-polarization, whereas SMEX02 observed~10% for V-polarization and~13% for H-polarization. Because of lower sensitivity to soil moisture and higher non-linear behavior of the radiative transfer model, these regions can be classified as heterogeneous environments. Due to the heterogeneous environment, there is a significant gradient observed in the sensitivity of land surface variables as they are upscaled linearly from 800 m to 12.8 km. However, as they are upscaled non-linearly, the gradient in the sensitivity measures decreases (Figures 5 and 6). Under heterogeneous environments, SMEX02 and SMAPVEX12 can further be classified as energy rich and water rich environments, respectively. SMEX02, which observed mean soil moisture~0.19 v/v and mean surface temperatures of~315 K, is classified as an energy rich environment. On the other hand, SMAPVEX12, which observed a relatively higher mean soil moisture~0.25 v/v and lower mean surface temperatures~290 K, is classified as water rich environment. For a heterogeneous water rich environment (SMAPVEX12), the sensitivity to land surface variables followed an exponential function with a support scale, and for a heterogeneous energy rich environment (SMEX02), the sensitivity to land surface variables followed a logarithmic function emphasizing the non-linearity/interactions of the region. Therefore, under heterogeneous environments, the upscaling method also determined the sensitivity to land surface variables.
In this study, a comprehensive framework is proposed based on the preliminary work presented by Neelam and Mohanty (2015) [26] (Figure 7), where the regions can primarily be classified as homogenous and heterogeneous environments based on density and heterogeneity in biomass. Each of these environments can further be classified as water and energy rich environments. The transition from homogenous to heterogeneous, energy rich to water rich, and vice versa, can occur through spatial parameters (extent and support), time (day, month, seasonality, etc.), and hydroclimates. For example, with a gradual increase in spatial scale (support, extent, and spacing) there will be a gradual change in land-cover/land-use, thereby changing heterogeneity in land surface variables captured with a varying spatial scale. Similarly, heterogeneity in land surface variables changes with temporal scales. For example, activities at monthly temporal scales such as seeding, irrigation, crop growth, harvesting, etc., and at seasonal scales such as spring, summer, fall, etc., effects heterogeneity in biomass, water, temperature, etc. Thus, the transitioning of the environments occurs through spatio-temporal scales, thereby changing the heterogeneity in the land surface variables and their sensitivities to the soil moisture retrieval algorithm. spatial parameters (extent and support), time (day, month, seasonality, etc.), and hydroclimates. For example, with a gradual increase in spatial scale (support, extent, and spacing) there will be a gradual change in land-cover/land-use, thereby changing heterogeneity in land surface variables captured with a varying spatial scale. Similarly, heterogeneity in land surface variables changes with temporal scales. For example, activities at monthly temporal scales such as seeding, irrigation, crop growth, harvesting, etc., and at seasonal scales such as spring, summer, fall, etc., effects heterogeneity in biomass, water, temperature, etc. Thus, the transitioning of the environments occurs through spatiotemporal scales, thereby changing the heterogeneity in the land surface variables and their sensitivities to the soil moisture retrieval algorithm. Figure 7. The conceptual model describes the classification of environments as homogenous and heterogeneous environments based on density and heterogeneity in biomass (green). Each of these environments can further be classified as water (blue) and energy (red) rich environments. Top-Left quadrant represents a heterogeneous water rich environment, where brightness temperature is most sensitive to vegetation and soil moisture with higher order interactions of ~5-10%; Top-Right quadrant represents a heterogeneous energy rich environment, where brightness temperature is most sensitive to vegetation and temperature with higher order interactions of >10-15%; Bottom-Left quadrant represents a homogeneous water rich environment, where brightness temperature is most sensitive to soil moisture with no or very low higher order interactions; Bottom-Right quadrant represents a homogeneous energy rich environment, where brightness temperature is most sensitive Figure 7. The conceptual model describes the classification of environments as homogenous and heterogeneous environments based on density and heterogeneity in biomass (green). Each of these environments can further be classified as water (blue) and energy (red) rich environments. Top-Left quadrant represents a heterogeneous water rich environment, where brightness temperature is most sensitive to vegetation and soil moisture with higher order interactions of~5-10%; Top-Right quadrant represents a heterogeneous energy rich environment, where brightness temperature is most sensitive to vegetation and temperature with higher order interactions of >10-15%; Bottom-Left quadrant represents a homogeneous water rich environment, where brightness temperature is most sensitive to soil moisture with no or very low higher order interactions; Bottom-Right quadrant represents a homogeneous energy rich environment, where brightness temperature is most sensitive to soil moisture and temperature with higher order interactions of <5-10%. The transition from homogenous to heterogeneous, energy rich to water rich, and vice versa, can occur through spatio-temporal scales (spatial: extent and support; time: day, month, seasonality, climate change, etc.) following land-use/land-cover change and events such as precipitation, evapotranspiration, etc. The lasting transition from homogeneous to heterogeneous and vice versa can also occur at the climate change temporal scales.

Conclusions
A comprehensive sensitivity analysis to study radiative transfer model using spatial maps is conducted for various field campaigns in multiple hydroclimates. The primary contribution of this work is to demonstrate how the sensitivity to the spatial maps of land surface variables change under various hydroclimates (Arizona, Oklahoma, Iowa, and Winnipeg) and evolving scales (0.8 km, 1.6 km, 3.2 km, 6.4 km, and 12.8 km) for a given extent. The impact of non-linear upscaling is illustrated by the comparison between linear and inverse distance weighted methods. The analysis resulted in environment-specific most sensitive variables: SM in homogeneous, and VWC and B in heterogeneous environments. Though the magnitude of sensitivity varied temporally, the ranking among the variables did not change within the study period. In the case of homogeneous landscapes or largely undisturbed natural environments, the sensitivity of T B to variables remained nearly uniform or followed an exponential function with support scales. Under higher heterogeneity, the T B sensitivity to vegetation and roughness followed a logarithmic function with an increasing support scale. Surface temperature always followed an exponential function under all conditions. The upscaling methods play a significant role in representing non-linear effects of heterogeneous landscapes, as processes and variables change their significance with scale. It is recommended that the most sensitive variables to T B (H-or V-) should be upscaled using non-linear (e.g., IDW) methods to preserve spatial heterogeneity, and this is especially important under heterogeneous environments. The study emphasized that the observed heterogeneity and upscaling method will determine the sensitivity to land surface variables. This study is particularly relevant in order to establish the significant variables that can be used for downscaling and upscaling T B algorithms to various scales and under heterogeneous landscapes.
The analysis assumed independence among the input variables, which might not be true under natural environments. However, there are very limited studies in the past that have established a correlation among land surface variables. Also, the study assumes the range of variability in the vegetation parameters and scattering albedo to be constant at all scales. Since there are no means to ascertain how those vary with scale, this assumption may introduce uncertainties. The future scope of this work can include analyzing the sensitivity of correlated land surface variables and its development with hydroclimates and scales.

Acknowledgments:
We acknowledge the support of the Schlumberger Faculty for Future Fellowship and thank all the students, scientists, and volunteers who participated in the various soil moisture campaigns and collected valuable data. The data for Oklahoma, Arizona, Iowa and Winnipeg, Canada used in this study can be accessed from the links provided below http://nsidc.org/data/amsr_validation/soil_moisture/index.html.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A
Consider a model defined as Y = f (X), where Y is the output, X = (X 1 , . . . X k ) are k independent inputs, and f is the model function. If all the factors X are allowed to vary over their range of uncertainty, then the corresponding uncertainty of the model output Y is defined by its unconditional variance V(Y). The input factors are ranked based on the amount of variance removed from the output when we learn the true value of a given input factor X i . That is, V X ∼i (Y X i = x * i ) would be the resulting conditional variance as it is fixed to its true value x * i taken over all X ∼i (all factors but X i ). However, the true value x * i is unknown for each X i . Hence, an average of V X ∼i (Y X i = x * i ) measured over all possible values x * i will remove the dependence of x * i , thereby resulting in E X i V X ∼i (Y|X i ) . Using the following property (Mood et al., 1950) [80], Thus, a small E X i V X ∼i (Y|X i ) , or a large V X i E X ∼i (Y|X i ) , will imply that X i is an important factor.
The conditional variance V X i E X ∼i (Y|X i ) is the first-order (e.g., additive) effect of X i on Y and the associated first-order sensitivity index of X i on Y measure, which is written as S i = V X i [EX ∼i (Y|X i )]

V(Y)
, while E X i V X ∼i (Y|X i ) is customarily called the residual. As such, S i is a number always between 0 and 1. Sobol [48] introduced the decomposition of model function f into summands of increasing dimensionality: In Equation (2), each individual term is a square integrable over the domain of existence and is mutually orthogonal [50]. These functions f i 1, ,...i k, are associated with partial variances through: and so on for higher order terms. Assuming all inputs to be independent, these terms are linked as: V(Y) = i V i + i j>i V ij + · · · + V 12...k . Dividing both sides of the equation by unconditional variance V(Y), we obtain; i S i + i j>i S ij + · · · + S 12...k = 1 (A5) where First Order Sensitivity Index S i = ; the amount of variance of Y explained by the interaction of the factors X i and X j (i.e., sensitivity to X i and X j not expressed in individual X i and X j ); Total Sensitivity Index S Ti = S i + i<j S ij + i< j<l S ijl + . . . S 1,2,...k ; this accounts for all the contributions to the output variation due to factor X i (i.e., first-order index plus all its interactions). S Ti can also be defined by decomposing the output variance V(Y), in terms of main effect and residual, conditioning with respect to all factors but one, i.e., X ∼i . The unconditional variance can be rewritten as: V(Y) − V X ∼i E X i (Y|X ∼i ) = E X ∼i V X i (Y|X ∼i ) , and dividing both sides by unconditional variance results in: . In other words, S Ti index is defined as a summation of the main, second, and higher order effects, which involves the evaluation over a full range of parameter space. For further details on the computation of sensitivity indices, please refer to Saltelli (2002).