Improving Estimates of Soil Salt Content by Using Two-Date Image Spectral Changes in Yinbei, China

: Soil salt content (SSC) is normally featured with obvious spatiotemporal variations in arid and semi-arid regions. Space factors such as elevation, temperature, and spatial locations are usually used as input variables for a model to estimate the SSC. However, whether temporal patterns of salt-affected soils (identiﬁed as temporal spectral patterns) can indicate the SSC level and be applied as a covariate in a model to estimate the SSC remains unclear. Hence, temporal changes in soil spectral patterns need to be characterized and explored as to their use as an input variable to improve SSC estimates. In this study, a total of 54 ﬁeld samples and a time-series of Sentinel-2 multispectral images taken at monthly intervals (from October 2017 to April 2018) were collected in the Yinbei area of western China. Then, two-date satellite images were used to quantify signiﬁcant spectral changes over time using spectral change vector analysis, and four two-date-based index methods were used to characterize soil spectral changes. Lastly, the optimal two-date-based spectral indices and multispectral bands were used as input variables to build the estimation models using a random forest algorithm. Results showed that the two-date-based spectral index could be applied as an input variable to improve the accuracy of SSC estimation at a regional scale. Temporal changes in salt-induced spectral patterns can be indicated by the band difference in the wavelength range from 400 nm to 900 nm. Three two-date-based indices designated as D 28a (i.e., the band difference between band 2 from an image acquired in April 2018 and band 8a from an image acquired in December 2017), D 22 , and D 28 were the optimal parameters for characterizing salt-induced spectral changes, which were dominated by the total brightness, chloride, and sulfate accumulation of the soils. The model did not yield satisfactory estimation results (RPD = 1.49) when multispectral bands were used as the input variables. Multispectral bands coupled with two two-date-based indices ( D 22 and D 28a ) used as the input variables produced the best estimation result (R 2 = 0.92, RPD = 3.27). Incorporating multispectral bands and two-date-based indices into the random forest model provides a remotely-sensed strategy that effectively supports the monitoring of soil salt content.


Introduction
Soil salinization is an important land degradation and desertification process that usually occurs in arid and semi-arid environments and poses a threat to the ecological balance. Soils containing excessive salt adversely impact the material exchange in green plant cells, and reduce grain yield and agricultural productivity [1,2]. Soil salinization is a worldwide issue resulting from anthropogenic activities (e.g., excessive agricultural practices and flood irrigation) and global climate warming [3,4]. Soil productivity for 40 to 45% of the global land surface area has decreased due to different levels of salinity on salt-affected soil [5].
Prioritization should be given to solving this difficult environmental issue by early recognition, timely detection, and accurate processing of data for monitoring total soil salt content (SSC) with little spatial variability [6]. Remote sensing has the advantages of wide-area coverage and the ability to use upgraded available soil datasets. Remote sensing has recently become an effective method for identifying and monitoring SSC levels in salt-affected soils [7]. Hyperspectral satellite data has been shown to be a powerful tool for SSC monitoring as result of refined spectral recognition for soil attributes, but it has had limited use for large-scale, high-frequency SSC detection because the acquisition costs of hyperspectral data are expensive, and because available data do not widely cover the world [8]. Multispectral satellite imagery is also one of the important remote sensing data sources, due to its temporal regularity, the large quantity of image sources, and low acquisition costs, and, therefore, could be used for SSC estimation and salinity monitoring [9][10][11][12][13][14].
SSC can be estimated from satellite data by directly using the reflectivity in multispectral data, or by checking indirect links between environmental elements (e.g., climatic attributes, organisms, and topography parameters) and saline soil based on the scorpan model theory [15,16]. In the first case, the regional scale spatial distribution of salt status can be estimated by the reflectivity of multispectral bands and simultaneously collected field data using multivariate regression methods [17][18][19][20]. In this method, the limited number of multispectral bands and broad wavelengths have not provided sufficient detection of SSC features, and have restricted the improvement in estimation accuracy. Furthermore, some scholars [21][22][23][24][25] have argued that the regression models using multispectral data-derived indices and field-collected physicochemical data can improve the mapping of SSC levels at the regional scale. The result is that the application of multispectral data-derived indices (e.g., soil salinity indices and various ecological indices) is limited by the factors of climate elements, atmospheric conditions, land cover, and the spectral resolution of the different satellite platforms. In the second case, the accuracy of estimation modelling has increased and produced promising results for monitoring SSC status at a regional scale because input covariates have included additional environmental factors, such as topography, geology, temperature, and moisture [1,[26][27][28][29][30][31]. This approach has great advantages in an area where the landscape has obvious variability. However, estimation accuracy would be decreased in large plain areas with a homogeneous environment.
As a matter of fact, apart from considering direct multispectral reflectance data and the indirect environmental influencing factors in the SSC model, the temporal spectral patterns are also a fundamental property of saline soils and are not sensitive to landscape patterns [16,32,33]. Specifically, salt accumulation in the soil is represented by a dense blocky formation or crust that is composed of white granular crystals on the surface layer. The combined crystal structure can increase the soil surface reflectance and produce significant spectral changes in the time of salt accumulation, such as changes in reflectance intensity and changes in the concave and convex pattern of spectral curves. These created spectral data show a close relationship with SSC level. Additionally, satellite images provide not only the geo-locations and spectral reflectance of the ground targets but also are able to provide temporal spectral information to illustrate the status of objects in study areas [34].
However, with regard to the temporal spectral information of salt-affected soils, whether it can indicate the SSC level, and whether it can be applied as an input variable in the model to estimate the SSC remains unclear. Hence, characterizing the temporal changes in soil spectral patterns and taking this as the input variable for improving estimates of SSC needs to be explored.
Characterizing the temporal changes in soil spectral patterns and modelling strategy selection are critical for improving SSC estimation. Two-date satellite images with significant spectral changes can contribute to providing temporal changes in spectral information that has a strong explanation capability regarding salt enrichment in the entire periods of salt accumulation. A two-date-based spectral index derived through the mathematical transformation of multidate satellite bands has been an effective tool for characterizing dynamic spectral information [35]. The two-date-based spectral index, including the band difference index, the ratio index, the normalized difference index, and the tangent of the angle formed by three bands, was constructed using paired multispectral satellite images to quantify the temporal changes in soil spectral information in this study. The difference index (D) and normalized difference index (ND) can be used to interpret spectral difference among two-date satellite bands because of the spectral intensity transitions of soil. Shape changes of spectral curves can be accessed by the tangent band index (tan). The ratio index (R) constructed by two-date satellite bands has the advantage of spectral feature integration and decreasing environmental impact in extracting soil spectral information. The generated spectral parameters reflect the different aspects of the physical and spectral change information of soil and can assist in the rapid recognition of salt-affected regions. However, an index that is powerful and strongly explains SSC levels still needs to be identified [21,[36][37][38]. An ensemble machine learning approach, random forest algorithm (RF), has been widely used in solving complicated regression and classification problems, and is becoming popular for optimal spectral parameter identification. RF can assess the global importance scores of each variable via the increasing class of regression trees [39], and the selected variables can give more comprehensive details regarding SSC and provide perspectives and understanding of the mechanisms involved with spectral parameters. In the regression process, the selected spectral indices used as input variables have reduced model complexity and enable a more effective determination of SSC status [40].
In the Yinbei region of western China, bare soils are permanently exposed to the remote sensors during the dry season, which provides the optimal time window for multispectral satellite estimation of SSC. In this severe drought period, agricultural management (e.g., fertilizer application, irrigation, and crop cultivation) is temporarily suspended; biophysical activities (e.g., humus transformation and soil texture variation) also slow down. In contrast, soil salt actively moves to the surface and accumulates continuously under dry conditions. Based on the above premises, we assumed that SSC level could be indicated by the temporal changes in salt-affected soil spectral reflectance. Hence, this study's objectives were to (1) estimate SSC using multispectral satellite images at a regional scale; (2) analyze spectral changes of two-date satellite images and build two-date-based spectral indices to characterize temporal changes in soil spectral reflectance; (3) evaluate the performance of the estimation model incorporating two-date-based spectral indices.

Materials and Methods
A flow chart of the analysis steps used in this study is shown in Figure 1. Acquisition and pre-processing of image data along with sampling and laboratory analysis were initially carried out. Secondly, spectral change analysis of multidate satellite images, construction of two-date-based indices, and determination of sensitive spectral parameters were completed. Lastly, RF models were determined and mapping results were obtained.

Study Area
The Yinbei region ( (Figure 2a). An alluvial plain is located in the center of the study area, and the Helan Mountain Range is located in the western part of the Yinbei region. The Yellow River on the eastern boundary of the region provides a large amount of irrigation water for agricultural practices in the central plain. The study area has a temperate continental climate with a mean temperature of 9.69 • C, a mean annual precipitation of 187 mm, and an average annual evaporation of 1774 mm. The elevation ranges between 1081 m and 3420 m above sea level. The irrigation-silted soil and the alluvial soil are the major soil type in the Yinbei area. The common geological structures are faulted basins with aquifers characterized by loose rock pore water. The Yinbei area is a typical arid region of western China that suffers from soil salinization due Remote Sens. 2021, 13, 4165 4 of 19 to environmental elements including high evaporation rate, high groundwater table, and poor drainage. Irrigation infiltration accounts for 80% of the total groundwater recharge, the proportion of evaporated water is over 76% of the total groundwater discharge, and the rest of the water comes from groundwater runoff [41]. These factors have intensified soil degradation and produced an extremely fragile ecological environment [42].
Remote Sens. 2021, 13, x 4 of 19 Figure 1. Flowchart of the process used for estimating soil salt content in the Yinbei region, China. RF: random forest; t0 and t1 represent optimal two-date satellite images; D: band difference; R: band ratio; ND: band normalized difference; tan: the tangent of the angle formed by three bands.

Study Area
The Yinbei region (38°35′56.4″ to 39°22′55.2″N and 105°58′12.0″ to 106°55′8.4″E) is located in the Northeast Ningxia Hui Autonomous Region of China, covers an area of 1000 km 2 , and has a population of 1.22 million (Figure 2a). An alluvial plain is located in the center of the study area, and the Helan Mountain Range is located in the western part of the Yinbei region. The Yellow River on the eastern boundary of the region provides a large amount of irrigation water for agricultural practices in the central plain. The study area has a temperate continental climate with a mean temperature of 9.69 °C, a mean annual precipitation of 187 mm, and an average annual evaporation of 1774 mm. The elevation ranges between 1081 m and 3420 m above sea level. The irrigation-silted soil and the alluvial soil are the major soil type in the Yinbei area. The common geological structures are faulted basins with aquifers characterized by loose rock pore water. The Yinbei area is a typical arid region of western China that suffers from soil salinization due to environmental elements including high evaporation rate, high groundwater table, and poor drainage. Irrigation infiltration accounts for 80% of the total groundwater recharge, the proportion of evaporated water is over 76% of the total groundwater discharge, and the rest of the water comes from groundwater runoff [41]. These factors have intensified soil degradation and produced an extremely fragile ecological environment [42].

Soil Sampling and Chemical Analysis
A total of 54 soil samples were collected from nine geographical areas (Figure 2b). Six samples were obtained within a 2 km by 1.5 km area at each site ( Figure 2c). If soil could not be collected at the designated sites due to being covered by river and road, a close location was substituted for sample collection. All sampling sites had their actual locations recorded with a handheld global positioning system. Due to the active processes of salt accumulation, most soil salt accumulated on the surface layer of soil (0-5 cm); thereby, the representative samples were collected based on this depth. All samples were collected in April 2018 and then sealed in polyethylene bags for later laboratory analysis. In the laboratory, all non-soil materials were removed from the samples, including stones and weeds, and then the samples were dried and passed through a 1 mm sieve under normal conditions (25 • C). Total SSC was measured using the weighted method in which soil and water were mixed at a ratio of 1:5, the soluble salt in the soil was extracted by agitation of the solution and filtration, the filtrate extraction was evaporated until dry, and then the soil organic matter was removed using sodium hydroxide [27]. Then, the SSC at each sampling location was determined.

Soil Sampling and Chemical Analysis
A total of 54 soil samples were collected from nine geographical areas ( Figure 2b). Six samples were obtained within a 2 km by 1.5 km area at each site ( Figure 2c). If soil could not be collected at the designated sites due to being covered by river and road, a close location was substituted for sample collection. All sampling sites had their actual locations recorded with a handheld global positioning system. Due to the active processes of salt accumulation, most soil salt accumulated on the surface layer of soil (0-5 cm); thereby, the representative samples were collected based on this depth. All samples were collected in April 2018 and then sealed in polyethylene bags for later laboratory analysis. In the laboratory, all non-soil materials were removed from the samples, including stones and weeds, and then the samples were dried and passed through a 1 mm sieve under normal conditions (25 °C). Total SSC was measured using the weighted method in which soil and water were mixed at a ratio of 1:5, the soluble salt in the soil was extracted by agitation of the solution and filtration, the filtrate extraction was evaporated until dry, and then the soil organic matter was removed using sodium hydroxide [27]. Then, the SSC at each sampling location was determined.

Time-Series Image Collection and Pre-Processing
The Sentinel-2 MultiSpectral Instrument (MSI) has a 10 m spatial resolution in four visible near-infrared spectroscopy bands, a 20 m spatial resolution in six red-edge and shortwave infrared (SWIR) bands, and a 60 m resolution for the coastal aerosol, water vapor, and cirrus band, the wavelength characteristics of bands are shown in Table 1 [43]. The revisiting frequency of Sentinel-2 MSI was up to 5 days. Monthly satellite images during the dry season (October 2017 to April 2018), including 21 cloudless (≤10%) Level-1 C sub-images, were acquired at no cost from https://glovis.usgs.gov/ (accesssed on 23 November 2020). Three MSI sub-images within a month were mosaicked into one representative image. Atmospheric correction was performed with Sentinel Application Platform software (SNAP, version 6.0) equipped with the Sen2cor plugin (version 2.8) that can transform MSI-L1C products with top of atmosphere reflectance to the MSI-L2A product with bottom of atmosphere reflectance values [44]. Band 1, band 9, and band 10 at 60 m resolution were used for aerosol calibration, water vapor correction, and atmospheric feature detection, respectively, and, thus, were not included in the subsequent

Time-Series Image Collection and Pre-Processing
The Sentinel-2 MultiSpectral Instrument (MSI) has a 10 m spatial resolution in four visible near-infrared spectroscopy bands, a 20 m spatial resolution in six red-edge and shortwave infrared (SWIR) bands, and a 60 m resolution for the coastal aerosol, water vapor, and cirrus band, the wavelength characteristics of bands are shown in Table 1 [43]. The revisiting frequency of Sentinel-2 MSI was up to 5 days. Monthly satellite images during the dry season (October 2017 to April 2018), including 21 cloudless (≤10%) Level-1 C sub-images, were acquired at no cost from https://glovis.usgs.gov/ (accessed on 23 November 2020). Three MSI sub-images within a month were mosaicked into one representative image. Atmospheric correction was performed with Sentinel Application Platform software (SNAP, version 6.0) equipped with the Sen2cor plugin (version 2.8) that can transform MSI-L1C products with top of atmosphere reflectance to the MSI-L2A product with bottom of atmosphere reflectance values [44]. Band 1, band 9, and band 10 at 60 m resolution were used for aerosol calibration, water vapor correction, and atmospheric feature detection, respectively, and, thus, were not included in the subsequent computation [11] The multispectral bands were all resampled at 10 m spatial resolution to reduce the adverse effects of mixed pixels. Additionally, the Yinbei region currently experiences periods of severe drought with low precipitation, and, therefore, the negative effect of soil moisture content on soil spectral reflectance is minimized.

Modelling Strategies
The remote sensing estimation model used the different spectral parameters as the input variables, as represented in Table 2. Selecting sensitive spectral bands (i.e., the reflectance reflectivity of the multispectral band is an important covariate) and optimal two-date-based spectral indices were vital parts of this modelling strategy. In the first, the sensitive spectral bands were characterized with the high Variable Importance Scores (VIS), that were identified by the random forest algorithm; the two-date-based index was then constructed. The VIS are also a key parameter for highlighting the optimal two-date-based index. Lastly, the remote sensing estimation models were established by the regression relationships between the input variables and measured SSC data. Table 2. List of the input variables and its sources for the random forest estimation model.

Number
Input Variables Data Sources Increasing SSC increases soil spectral reflectance resulting from the high reflectivity of white salt granules. Two-date satellite images indicating obvious spectral changes are useful for revealing the salt increment on the surface soil layer in the entire period of salt accumulation [45,46]. Spectral change vector analysis is an effective method for evaluating the direction and magnitude of soil spectral reflectance changes that result from the changes in temporal salinization information, and facilitate the identification of the two-date images with significant spectral change information [47,48]. The computation process is given in Equation (1).
where g i and h j represent the spectral vector of the image pixels of the temporal ith image and the temporal jth image, respectively, ∆G is the change intensity value. The calculation was conducted using MATLAB 2017b (MathWorks, Natick, MA, USA). Soil spectral changes for the multidate satellite images included six two-date images (October 2017 and April 2018, November 2017 and April 2018, January 2018 and April 2018, February 2018 and April 2018, and March 2018 and April 2018), and were analyzed using spectral change vector analysis. The soil's spectral reflectance changes in six two-date images can be judged by the change intensity value ∆G , with high values suggesting significant changes in soil spectral reflectance. The simultaneous temporal images captured in April 2018 were applied in concert with the field samples. After harvest, bare soils covered the Yinbei area beginning in October. Therefore, acquisition of the time-series images started in October 2017. Due to the availability of cloudless images for the Yinbei area, detection of soil spectral reflectance was conducted monthly (every 30 days). Additionally, the minimal period between images for capturing spectral reflectance changes due to soil salinization depended on the local area's environmental condition (e.g., groundwater level, precipitation, and evaporation rate). In the study area, salt accumulation is a continuous process in long-term drought conditions; monthly satellite images can meet the requirement for detecting the temporal changes in soil spectral patterns.

Two-Date-Based Spectral Index Construction
The two-date-based spectral index is an effective tool that can help in acquiring temporal spectral patterns from the saline soil system [34,35]. As with the determination of the SSC sensitive bands from the simultaneous image and the optimal two-date images (t 0 and t 1 ), construction of the two-date-based index was based on the D, R, ND, and tan approaches. The equations for index construction are provided in Table 3. Figure 3 presents a case of two-date-based index establishment. The mathematical transformations of satellite bands were calculated with ENVI 5.3 IDL software (Exelis Visual Information Solutions, Boulder, CO, USA). Table 3. Two-date-based spectral index construction using the two-date images (t 0 and t 1 ) in the Yinbei region, China.

Index Abbreviation
Calculation Image Equation Note: B i is the ith sensitive band in the multispectral image (t 0 ), B j is the jth sensitive band in the image (t 1 ), D ij represents the difference between B i and B j , R ij indicates the ratio of B i and B j , ND ij is the normalized difference index of band i and band j; tan imj is the tangent of the angle calculated by three bands (ith band, mth band, and jth band), where the mth band comes from (t 0 ) and the other bands are derived from (t 1 ). k im is computed by the slope of the linking line between ith and mth band. k mj is also calculated by the slope of the linking line between mth and jth band.
temporal spectral patterns from the saline soil system [34,35]. As with the determination of the SSC sensitive bands from the simultaneous image and the optimal two-date images (t0 and t1), construction of the two-date-based index was based on the D, R, ND, and tan approaches. The equations for index construction are provided in Table 3. Figure 3 presents a case of two-date-based index establishment. The mathematical transformations of satellite bands were calculated with ENVI 5.3 IDL software (Exelis Visual Information Solutions, Boulder, CO, USA). Table 3. Two-date-based spectral index construction using the two-date images (t0 and t1) in the Yinbei region, China.

Index Abbreviation Calculation Image Equation
Note: Bi is the ith sensitive band in the multispectral image (t0), Bj is the jth sensitive band in the image (t1), Dij represents the difference between Bi and Bj, Rij indicates the ratio of Bi and Bj, NDij is the normalized difference index of band i and band j; tanimj is the tangent of the angle calculated by three bands (ith band, mth band, and jth band), where the mth band comes from (t0) and the other bands are derived from (t1). kim is computed by the slope of the linking line between ith and mth band. kmj is also calculated by the slope of the linking line between mth and jth band.  R 1212 is the ratio in B12 between t 0 and t 1 ; ND 28a is the normalized difference index of band 2 and band 8a. tan 288a represents the tangent of the angle of three bands (B2, B8, and B8a) of images t 0 and t 1 , where B8 comes from image t 0 and B8a and B2 come from image t 1 .

Random Forest Algorithm
In the random forest (RF) computation, the kth sub-datasets are randomly sampled from the original dataset, the β features are extracted from the m features of the sample, and the sub-dataset is modelled and optimized by the rest of the dataset [49]. Variable Importance Scores (VIS) depended on the proportion of the increment in the mean square error. In other words, model accuracy was assessed by the out-of-bag (OOB) error. A high VIS value means that the close relationship between spectral parameters and SSC, as well as the spectral parameters, have an obvious influence on the estimation result. VIS is expressed as follows where t is the tree in the forest, OOB t is the linked portion not involved in bootstrap procedures for t establishment, errOOB t shows the error of the tth tree linked with the OOB t sample, and errO OB j t is the disturbed samplings of the X j permuted value. Input parameter selection and model estimation used the RF algorithm mentioned above and can be classified into three steps. In the first step, the collected spectral parameters (multispectral bands or spectral indices) were assessed for their ability to interpret SSC using the RF algorithm. Then, the same procedures were re-implemented to compute the VIS after discarding the 10% of parameters having the lowest VIS values until the number of spectral parameters in the remaining subset was less than 10 [50]. In order to ensure that there were a sufficient number of input variables and to save training time, a VIS threshold (10%) was set to identify the optimal spectral parameters [34]. Lastly, optimal spectral parameters were selected and participated in the subsequent calculations for model training and validation. The size of the input data (mtry) and the quantity of trees (ntree) were set to 7 and 1000, respectively, in the calculation. The VIS value and RF algorithm were processed via the R statistical computing package (R version 4.0.2, R Development Core Team, 2020). More details can be found in [51] and [49]. In order to test the applicability and the stability of the estimation model, the soil salt contents of the Yinbei region were estimated and mapped using the optimal constructed model.

Statistical Analysis
Based on the equations in Table 3, 55 two-date-based indices were generated based on the transformation of sensitive bands among two-date images, including 15 D indices, 15 R indices, 15 ND indices, and ten tan indices.
The coefficient of determination (R 2 ) and the ratio of performance to deviation (RPD) were the statistical parameters for evaluating model accuracy and stability. The R 2 value is the proportion of the target variable variation that is explained by the model, and the RPD value is the ratio of the standard deviation to the root mean square error. In general, the best model is characterized by high R 2 and RPD values [52]. In the calibration and validation processes, a leave-one-out cross-validation was applied to assess model accuracy and stability, and the total dataset at each sampling involved the calibration and validation procession. Specifically, one sample was part of the calibration set after validation competition. According to the classification criterion from [53], RPD > 2 suggests that the model is more reliable with a higher predictive power, RPD values between 1.4 and 2 indicate that the model can generate moderately accurate outcomes, and RPD < 1.4 suggests that the model does not have adequate predictive capabilities.

Characteristics of Soil Properties and Reflectance Spectra
The average value of SSC in the surface layer (0-5 cm) was 44.59 g kg −1 , with values of 2.92 g kg −1 to 290.86 g kg −1 . The soil dataset had significant spatial variability in total salt contents, with the coefficient of variation for total salt content reaching 1.34, indicating a high variability level and the existence of representative values from the different soil type. The mean pH was 8.56, suggesting that the Yinbei area was moderately alkaline affected. The coefficient of skewness for SSC was 2.97. The distribution of elemental contents should possess a normal distribution, and a positive skew suggests that the outlier values arose from an exterior disturbance. According to the statistical analysis of SSC status, the sampling sites exhibited a high variability level of SSC, and a small portion of the hotspot area containing a high salt content represented an active process of soil salt accumulation.
The spectral reflectance visually showed the variation between samples with different levels of SSC (Figure 4). Reflectance intensity of the soil from sampled locations increased as the salt content rose from 2.92 g kg −1 to 290.86 g kg −1 in the MSI dataset. As shown in Figure 4, high salt content soils are more likely to exhibit a white salt crust on the topsoil and are characterized by higher spectral reflectance than soils with a lower soil salt content status. The correlation coefficients between satellite spectral data and measured SSC are displayed in Table 4. Reflectance intensity of all bands exhibited a positive correlation with SSC; band 2, band 8, band 8a, band 11, and band 12 were significantly correlated with SSC. Specifically, the reflectance of five bands changed significantly in soil with different salt-affected levels, showing higher reflectance for soils with SSC contents of 113.41 g kg −1 and 290.86 g kg −1 than for soil with SSC contents of 2.92 g kg −1 and 31.02 g kg −1 . Soil reflectance increased as soil salt content increased across the entire spectrum, and the high spectral reflectance areas were located in the near-infrared wavelengths. In contrast, low spectral reflectance was observed in the visible spectra. The spectral curves shown in Figure 4 showed rapidly increasing reflectance in the visible regions as salt content increased. In general, the clearly different reflectance intensities of the spectral curves in the visible and near-infrared regions were helpful for distinguishing soils having different salt contents.

Optimal Two-Date Satellite Images
The results of change vector analysis are shown in Figure 5 ( Figure 5a represents the conceptual graph of changes of spectral vectors, and Figure 5b is the change intention values of different two-date images). The temporal changes in the spectral reflectance from time-series images (October 2017 to March 2018) to the simultaneous temporal image (April 2018, t0) have obvious differences, as shown by the change intensity value (‖ΔG‖) in the range from 0.04 to 0.34. The change intensity value between the December vector and the April vector reached 0.34, indicating that the spectral reflectance change for saline soils was highest. The January vector and the April vector also had a high change intensity value (0.31) that captured the marked changes in temporal spectral information related to soil salinization. The other two-date images did not show significant changes in the spectral reflectance of the soil as the change intensity value was lower than 0.30. In order to extract temporal changes in soil spectral information in the whole period of salt accumulation, the largest change intensity value was regarded as a standard for selecting the optimal two-date images. Therefore, the ideal times for acquiring two-date images (t0 and t1) were determined to be in April 2018 and December 2017, respectively.  Changes in soil spectral reflectance showed a close relationship with SSC level. Dynamic hydrological processes drove salt accumulation in the topsoil of the Yinbei region, and created noticeable changes in the spectral reflectance of field soil. Four stages of the time of soil salt accumulation are illustrated in Figure 4. In the left panel, layer (a) shows a healthy status soil (low salt) after the crop is harvested and is characterized by grayness and lower reflectance. Following the salinization process, soil brightness increases (layer (b)). Layer (c) shows that the soil has been moderately affected by salt. Soil spectral reflectance at this stage began to represent salt mineral characteristics as a white salt crust begins to form. Lastly, salt signals significantly influence the soil spectral pattern of image pixels (layer (d)), showing the highest brightness and reflectivity among salt-affected soils. The spectral reflectance of the saline soil changes significantly in the time of active salinization processes (from layers (a) to (d)). Increasing soil surface brightness can enhance the soil's reflecting power of electromagnetic energy and alter the reflectance intensities of the satellite spectral band. Additionally, the reflectance values in the range from 950 nm to 1650 nm were dominated by the stretch bending of the C-O and COO bonds; OH-Al banding form has a close relationship with the clay minerals and exhibited absorption features at 2200 nm; the mineral gypsum (CaSO 4 ) is also an important component of saline soil, and spectral reflectance was observed at approximately 1730 nm. The white salt granules may be related with internal vibration in CaCO 3 and the generated absorption bands around the wavelength of 2350 nm [31]. Thus, salt accumulation exhibited changes in spectral reflectance (reflectance intensity, changes in the concave and convex pattern of spectral curves). These spectral data show a close relationship with SSC and can be utilized as an SSC feature to provide the basis for improving the SSC estimation model.

Optimal Two-Date Satellite Images
The results of change vector analysis are shown in Figure 5 (Figure 5a represents the conceptual graph of changes of spectral vectors, and Figure 5b is the change intention values of different two-date images). The temporal changes in the spectral reflectance from time-series images (October 2017 to March 2018) to the simultaneous temporal image (April 2018, t 0 ) have obvious differences, as shown by the change intensity value ( ∆G ) in the range from 0.04 to 0.34. The change intensity value between the December vector and the April vector reached 0.34, indicating that the spectral reflectance change for saline soils was highest. The January vector and the April vector also had a high change intensity value (0.31) that captured the marked changes in temporal spectral information related to soil salinization. The other two-date images did not show significant changes in the spectral reflectance of the soil as the change intensity value was lower than 0.30. In order to extract temporal changes in soil spectral information in the whole period of salt accumulation, the largest change intensity value was regarded as a standard for selecting the optimal two-date images. Therefore, the ideal times for acquiring two-date images (t 0 and t 1 ) were determined to be in April 2018 and December 2017, respectively.

Sensitive Spectral Parameters
The RF importance scores of 10 single bands are shown in Figure 6a. B11 (shortwave infrared band) showed significant importance for soil salt. B2, B8, B8a, and B12 also have high VIS values (>10 %) and were selected as sensitive bands. The variable importance assessment approach was also implemented to select the optimal indices to quantify the salt-induced spectral changes. The calculation results are shown in Figure 6b. Obvious differences in VIS for the two-date-based indices existed in the four transformation methods. The most important indices (VIS > 10%) came from the band difference indices, such as D22, D28a, and D28. The ratio index scores did not rank in the top 10 and were lower than most of the index transformed methods calculated by the band difference, the band normalized difference, and the bands tangent. The band difference index was the most important method for highlighting salt-induced spectral changes. Therefore, five sensitive bands (B2, B8, B8a, B11, and B12) were selected, and the three optimal two-date-based indices (D22, D28a, and D28) were employed to characterize the temporal spectral patterns of salt-affected soils.

Sensitive Spectral Parameters
The RF importance scores of 10 single bands are shown in Figure 6a. B11 (shortwave infrared band) showed significant importance for soil salt. B2, B8, B8a, and B12 also have high VIS values (>10 %) and were selected as sensitive bands. The variable importance assessment approach was also implemented to select the optimal indices to quantify the salt-induced spectral changes. The calculation results are shown in Figure 6b. Obvious differences in VIS for the two-date-based indices existed in the four transformation methods. The most important indices (VIS > 10%) came from the band difference indices, such as D 22 , D 28a , and D 28 . The ratio index scores did not rank in the top 10 and were lower than most of the index transformed methods calculated by the band difference, the band normalized difference, and the bands tangent. The band difference index was the most important method for highlighting salt-induced spectral changes. Therefore, five sensitive bands (B2, B8, B8a, B11, and B12) were selected, and the three optimal two-date-based indices (D 22 , D 28a , and D 28 ) were employed to characterize the temporal spectral patterns of salt-affected soils.
methods. The most important indices (VIS > 10%) came from the band difference indices, such as D22, D28a, and D28. The ratio index scores did not rank in the top 10 and were lower than most of the index transformed methods calculated by the band difference, the band normalized difference, and the bands tangent. The band difference index was the most important method for highlighting salt-induced spectral changes. Therefore, five sensitive bands (B2, B8, B8a, B11, and B12) were selected, and the three optimal two-date-based indices (D22, D28a, and D28) were employed to characterize the temporal spectral patterns of salt-affected soils.

Performance Assessment of Prediction Models Built with Different Inputs
Results of the estimation models are shown in Figure 7. The models built with five sensitive bands (B2, B8, B8a, B11, and B12) and three optimal two-date-based indices (D 22 , D 28a , and D 28 ) had RPD values of 1.49 and 1.33, respectively, indicating that the models did not produce satisfactory estimation results. The RF model was constructed when five sensitive bands combined with an optimal two-date-based index (D 22 , D 28a , and D 28 ) as an input variable, as represented by the RPD values of 2.75, 2.52, and 2.41, respectively (Figure 7c-e). The optimal two-date-based indices used as a part of the input variable showed better accuracy than the five sensitive bands in the RF model. As shown in Figure 7f-i, different indices combinations were created to further analyze the capability of two-date-based indices to improve model accuracy. Figure 7c-h suggest that the model performance increased when the optimal two-date-based indices combination was used as an auxiliary input, as RPD values improved from the range of 2.41 to 2.75 to the range of 2.72 to 3.27. Five sensitive bands (B2, B8, B8a, B11, and B12) coupled with two optimal two-date-based indices (D 22 and D 28a ) in the RF algorithm produced satisfactory estimation results (R 2 = 0.92, RPD = 3.27) with the highest accuracy and stability.

Map of Soil Salt Contents
The model was applied to each pixel of the MSI dataset in the study area, and the mapping result is presented in Figure 8. The results indicated that strong salt-affected areas were located in the southwest Yinbei area and on the two sides of the Yellow River, accounting for the areas of highest evaporation and lowest depth to groundwater. The growing demand for agricultural irrigation promotes salt water transportation and results in a considerable accumulation of soluble salt in surface soil. Figure 8 shows that SSC values in most of the Yinbei area were lower than 31.40 g kg −1 (72.06%). The results similarly identified the aforementioned descriptive statistics of this region. Additionally, high salt soils (31.40 g kg −1 < SSC ≤ 302.40 g kg −1 ) were observed to be scattered throughout the entire area shown in Figure 8. This can be explained by the fact that soil microtopography is shaped by agricultural irrigation that enables farmland to hold large amounts of irrigation water and leads to rising groundwater levels. Irrigation coupled with high evaporation results in a high level of SSC in the dry season. two-date-based indices to improve model accuracy. Figure 7c-h suggest that the model performance increased when the optimal two-date-based indices combination was used as an auxiliary input, as RPD values improved from the range of 2.41 to 2.75 to the range of 2.72 to 3.27. Five sensitive bands (B2, B8, B8a, B11, and B12) coupled with two optimal two-date-based indices (D22 and D28a) in the RF algorithm produced satisfactory estimation results (R 2 = 0.92, RPD = 3.27) with the highest accuracy and stability.

Map of Soil Salt Contents
The model was applied to each pixel of the MSI dataset in the study area, and the mapping result is presented in Figure 8. The results indicated that strong salt-affected areas were located in the southwest Yinbei area and on the two sides of the Yellow River, accounting for the areas of highest evaporation and lowest depth to groundwater. The growing demand for agricultural irrigation promotes salt water transportation and results in a considerable accumulation of soluble salt in surface soil. Figure 8 shows that SSC values in most of the Yinbei area were lower than 31.40 g kg −1 (72.06%). The results similarly identified the aforementioned descriptive statistics of this region. Additionally, high salt soils (31.40 g kg −1 < SSC ≤ 302.40 g kg −1 ) were observed to be scattered throughout the entire area shown in Figure 8. This can be explained by the fact that soil microtopography is shaped by agricultural irrigation that enables farmland to hold large amounts of irrigation water and leads to rising groundwater levels. Irrigation coupled with high evaporation results in a high level of SSC in the dry season.

Characterizing Temporal Changes in Salt-Induced Spectral Information for Improving SSC Estimation
Choosing sensitive wavelength range was the key step for characterizing temporal changes in salt-induced spectral information. Several previous studies [11,18,54] have shown the close relationship between shortwave infrared bands and SSC, and regarded them as the key input parameters in SSC estimation. Nevertheless, the selected two-date-based indices did not involve the combination with B11. Bannari et al. [22], Jin et al. [55], and Shrestha [56] have suggested that band 11 is the critical indicator of soil salt. Areas with increased soil salinity levels normally exhibit increased coverage of the white salt crust, causing increased spectral reflectivity in almost every multispectral band (400-2500 nm) [11]. The different thicknesses of white salt crust in surface soils are mainly chloride and sulfate in the Yinbei area. Sodium chloride is normally observed as white granular crystals in the soil. In contrast, sodium sulfate belongs to the monoclinic crystal system and is observed as a dense blocky or shell-like structure on the soil surface [46,[57][58][59]. The result is that high brightness with increasing reflectivity is produced in almost every spectral band of the soil. At the same time, the salt accumulation process in surface soils, soil alkalization, soil sodicity, and whitening of the surface color will also occur [60,61]. Bai et al. [17], Day [62], Lopez-Granados et al. [63], and Stoner et al. [64] have proven the significant relationships between SSC and soil alkalinity, with salt accumulated in soils that resulted from high pH and decreased soil organic matter. These changes influence soil color and the behavior of spectral features in the wavelength range from 350 nm to 1200 nm, and can be regarded as critical indicators of salt in field soil. The spectral bands responding to the salt-induced spectral changes in this study have also been observed at the susceptible regions from the visible near-infrared range (400-900 nm), similar to our findings. Therefore, the dynamic spectral information associated with SSC were especially sensitive to the band difference in the visible near-infrared range (400-900 nm) for the two-date satellite images.
The optimal estimation model can be constructed by the two-date-based indices with high VIS values. As shown in Figure 6, the VIS value of D 22 was higher than the value of D 28a of 10.55%; the estimation errors were 8.36% lower than for the model using D 28a as an auxiliary input. Additionally, increasing the number of two-date-based indices used as input variables had a positive effect on the improvement in the model's estimation accuracy, although a little estimation error of two-date models as a result of information redundancy between two-date-based indices will also arise. As shown in Figure 7c to Figure 7h, the accuracy of the two-date images model with seven independent variable inputs (five sensitive bands coupled with two two-date-based indices) was greater than that of the model with six independent variable inputs (five sensitive bands coupled with one two-date-based index). Nevertheless, the best estimation model was achieved when the auxiliary inputs were two two-date-based indices (D 22 and D 28a ) rather than three twodate-based indices (D 22 , D 28 and D 28a ). Collinear information existed between D 28a and D 28 because B8 and B8a contained common wavelength ranges and had similar responses for soils. Therefore, estimation errors arose for the model built with three two-date-based indices (D 22 , D 28 and D 28a ).
The absence of a temporal spectral analysis of soil could restrict the performance improvement of an estimation model based on multispectral images. The results presented in this study provide a solution for analyzing the temporal changes in soil spectral patterns and increasing the explanatory ability of a multispectral remote sensing model for SSC, thereby promoting the ability to estimate SSC in an arid land.

Robustness of the Model Built with Two-Date-Based Indices
A satisfactory result was obtained using the two-date-based indices from the optimal two-date satellite images. The other two-date images during the bare soil period were employed to further examine the capability of the two-date image-based models. The optimal estimation results of the remaining five two-date image-based models are provided in Table 5. The two-date-based indices calculated from the remaining five two-date images also produced positive effects on the estimation accuracy of SSC modelling. The SSC model built with two-date images always showed better accuracy than the result derived from the single simultaneous image (Figure 7a). The two-date-based index's introduction to the estimation model significantly improved accuracy. For example, the RPD values for the twodate model (January 2018 and April 2018) increased from 1.49 to 2.72. However, there were performance differences in the five two-date models. Specifically, the significance of saltinduced spectral changes was different between the six two-date images. The explanation capability of two-date-based indices for SSC was different, affecting the performance of two-date models. The best result and lowest accuracy of the two-date model were produced by the two-date images with obvious soil spectral changes (December 2017 and April 2018) and insignificant spectral changes (October 2017 and April 2018), respectively. Additionally, the multi-temporal images (three date images or more) could provide extra spectral covariates (i.e., multidate-based spectral indices) and could result in high degrees of information redundancy among spectral indices. Therefore, the effects on SSC modelling performance should be explored in the future. Note: R 2 is the coefficient of determination; RPD is the ratio of performance to deviation.
The modelling strategies were also extended to the Operational Land Imager (OLI) dataset to verify the capability of the model built with the two-date-based indices. All bands in the OLI dataset exhibited a strong linear relationship with the MSI spectral reflectance on ground samples [54]. Two cloudless OLI images (December 2017 and April 2018) were acquired from https://glovis.usgs.gov/ (accessed on 25 November 2020) at no cost. Pre-processing procedures, including geometric correction, radiation correction, and atmospheric correction, were then performed [65]. The calculation results are shown in Table 5. RF using four sensitive bands (band 4, band 5, band 6, and band 7) as independent variables did not provide acceptable accuracy of SSC estimation (i.e., RPD value was lower than 1.40). The model's accuracy was increased when the two-date-based indices (D 55 ) were introduced as input variables, as the RPD value increased from 1.25 to 1.91. Results indicated that D 55 was the only index with VIS value (37.79%) higher than 10%, and D 55 was still the difference index combined by the bands in the near-infrared wavelengths of 400-900 nm. These findings verify the results from the MSI datasets that the saltinduced spectral change was more sensitive to different bands in the wavelength range of 400-900 nm among two-date images. The difference band index may be an effective method for quantifying the changes in temporal spectral patterns in the case of salt-affected soils. Analysis of the physical signals and temporal spectral patterns of objects in multispectral images boosts the SSC model's accuracy and stability, and the generated accuracy was approximately equal to the accuracy of the models built with hyperspectral data reported in previous studies [66][67][68].
In general, temporal spectral patterns are a basic property of saline soils [33]. Extracted spectral change information (identified as the two-date-based spectral indices) can be an indicator of the SSC level, and, as such, we introduced them as input variables that increased the explanatory ability of estimation models for SSC. This modelling strategy has the benefits of simple operation and decreasing field work associated with investigating environmental covariates of salinization. This modelling strategy was also more appropriate for application in arid areas to which travel is inconvenient and climate is harsh. The temporal spectral patterns of salt-affected soil can also be used as a proxy with a new index design for monitoring SSC status.

Uncertainty Analysis
Salt accumulation in the surface soil and the associated temporal changes in the spectral patterns may not be linearly related. Apart from the aforementioned temporal saltinduced factors that can significantly affect the spectral reflectance of soils, environmental factors such as vegetation and moisture will also produce a certain influence on reflection intensity and the shapes of the spectral curve [69,70]. However, this part of the impact is minimized in the dry season of the Yinbei region of western China. The spectral changes in the period of salt accumulation constitute a complex synthetic signal that consists of soil substance signals (e.g., soil organic matter and soil irons), soil texture signals, and environmental signals [71]. As such, an analysis of the mechanisms causing salt-induced spectral changes will require further strict laboratory experiments.
The cross-validation method for remote sensing estimation of SSC may be optimistically biased. Although a leave-one-out cross-validation is an acceptable option for quality measures of the soil map, the autocorrelation existed in the soil dataset that can lead to overestimated estimation accuracy [72]. Designing an unbiased validation method and a design-based sampling strategy may be ideal for quality measures of remote sensing estimation [73]. However, that is beyond the range of our current study, and deserves to be studied in the next.
Mixed pixels of satellite images and small size datasets of field samples were the important factors contributing to modelling uncertainty. The small size dataset of field samples may increase data distribution variation and improve regression complications. Apart from ensuring the representation of each sampling collection, the RF model is a fully non-parametric method and did not assume linear linking in input variables, thereby promoting estimation efficiency such that a large sample size is not necessary [74,75]. Mixed pixels in the imagery dataset are an unavoidable problem [76]. Including the nonsoil end-members leads to misidentification of spectral features and increased uncertainty. Soil sampling locations were designed to be at sites with pure soil pixels by assessing geographical consistency. Non-soil regions were excluded using the mask method to reduce the negative effect of outliers. The images mask is also implemented to decrease the interference of a foreign body with the spectrum, especially a few sandy soil areas on both sides of the Yellow River that had a spectral reflectance similar to that of saline soil. Therefore, ensuring that a soil's reflectance is not contaminated by other targets is important for obtaining accurate results, and the decomposition of mixing pixels is still a topic to be studied in the future.

Conclusions
Characterizing temporal changes in soil spectral patterns as an input variable for an RF model can improve the remote sensing estimation of SSC (0-5 cm) at a regional scale. Results showed that band difference in the visible near-infrared range (400-900 nm) among two-date satellite bands was an effective method for characterizing dynamic spectral reflectance of the salt-affected soils. The total brightness, chloride, and sulfate accumulation of soils dominated optimal spectral parameters for characterizing salt-induced spectral changes. The sensitive multispectral bands included B2, B8, B8a, B11, and B12. B11 had the highest variable importance score. Sensitive bands as model inputs did not produce satisfactory estimation results (RPD = 1.49). Estimation accuracy for the model was increased when the two-date-based index was included as an input variable. The model built with the combination of five sensitive bands (B2, B8, B8a, B11, and B12) and two two-date-based indices (D 22 and D 28a ) produced the best prediction result (R 2 = 0.92, RPD = 3.27). Accurate and timely identification of changes in SSC on arid and semi-arid lands by use of the method we have described in this study can help producers and policy makers improve the understanding of the current salt-affected status so that they can implement effective measures for land restoration and soil reclamation.

Data Availability Statement:
The data presented in this research can be requested from the corresponding author or first author.