Analysis on the Spatio-Temporal Changes of LST and Its Inﬂuencing Factors Based on VIC Model in the Arid Region from 1960 to 2017: An Example of the Ebinur Lake Watershed, Xinjiang, China

: LST (Land surface temperature) is an important indicator for monitoring dynamic changes in the earth’s resources and environment. However, the complexity of obtaining long-term, continuous LST data hinders the development of research on LST responses to meteorological factors or LUCC in areas where data is lacking. The objective of this research was to use the VIC-3L (Variable Inﬁltration Capacity) based on multi-source remote sensing data to simulate and explore spatio-temporal changes in the LST, to analyze the relationship between the LST and meteorological elements by using cross-wavelet transform (XWT) and wavelet coherence (WTC), the relationship between the LST and LUCC by using three-phase remote sensing images of LUCC. The following results were obtained. The annual average LST of the study area is increasing at a rate of 0.027 ◦ C per year. The annual average LST level is relatively high in the central and eastern regions. The average temperature has an important inﬂuence on LST, which is mainly reﬂected in the period scale of 1~4a in 1963–1972, 1980–1996, and 2004–2010. The sharp decline in open shrubs may have exacerbated the increase in LST in the study area. This study provides a scientiﬁc reference for studying LST in arid areas.


Introduction
Land surface temperature (LST) is a typical indicator of the energy flow in the land surface-atmosphere-biosphere interactions. The UN Intergovernmental Panel on Climate Change (IPCC) Fifth Assessment Report pointed out that the global average surface temperature increased by 0.89 • C (0.69~1.08 • C) from 1901 to 2012 [1]. The change of LST is an intuitive regional climate response to global climate change. On the one hand, the LST variations will directly affect the soil moisture content and the soil evaporation, thereby indirectly affecting regional vegetation growth, agricultural production, and the water cycle [2,3]. On the other hand, the LST is an important element for predicting drought, climate change [4][5][6]. Therefore, it is of greater research significance in agriculture, hydrology, ecology, environment, climate, and biogeochemistry [4,[7][8][9][10][11].
Research on changes in LST is more important in arid regions. One of the most prominent features of arid regions is the low atmospheric precipitation [12], which directly affects the soil moisture content and indirectly affects the vegetation conditions in these areas [13,14]. An increase in the LST can lead to reduce soil moisture and increase dust emission [15]. LST can affect wind erosion through changes in soil moisture content, as temperature increases, the evapotranspiration from the land surface increases, and the soil moisture content decreases [16]. As a result, the adhesion between the soil particles is reduced, and the soil particles are more easily separated from the soil, thereby increasing the concentration of the particles in the atmosphere [17]. Moreover, LST has shown a positive relationship with changes in the frequency of wind erosion events in Asia [18]. Therefore, studying the temporal and spatial changes of the LST and its driving factors is helpful to better understand the evolution of drought and water resources and provide theoretical support for sustainable development.
Under the global warming background, experts and scholars in different fields are paying more and more attention to changes in LST. The evolution of the global LST trend in the past century is diagnosed using the spatial-temporally multidimensional ensemble empirical mode decomposition method, the fastest warming in recent decades (>0.4 K per decade) occurred in northern mid-latitudes [19]. It also showed a particularly enhanced trend of LST in semi-arid regions [20]. However, the existing studies on LST dynamics mainly focused on the depiction and quantitative temporal comparison of LST spatial distribution or LST time series analysis over time [21]. The relationship between LST and other climate elements is less involved. Moreover, there are still some difficulties in the long-term large-scale continuous simulation of LST in data-scarce areas.
Exploring the mechanism of interaction between LUCC and local climate has become a popular research topic [22][23][24]. In the past few years, many studies based on the impact of LUCC on LST have focused on the urban heat island (UHI) [25][26][27]. The increase of the impermeable layer increases the absorption of solar radiation to a certain extent, thereby increasing the LST [28]. However, each land cover type possesses unique qualities in terms of energy radiation and absorption, in addition, LUCC also alters greenhouse gas emissions, such as CO 2 , in the atmosphere [29]. Although the urbanization process in some areas is not fast, the effect of LUCC on LST is worth studying, especially for arid areas.
The Land surface model (LSM) is a key tool to predict the terrestrial hydrosphere and its atmospheric coupling and to understand this in-depth. Distributed LSM predicts hydrological states and fluxes, such as LST, at each grid cell [30]. The LSM based on multisource remote sensing data can carry out large-scale simulations of LST to compensate for the discontinuity in time and space caused by only using remote sensing images and station data. The VIC (Variable Infiltration Capacity) model is one of the core application models of the GLDAS and NLDAS [31]. Therefore, it is effective to use VIC to carry out large-scale and long-term continuous LST simulations in scarce data areas.
The main objectives of this study are (1) to simulate daily LST from 1960 to 2017 by using multi-source data and to verify the accuracy of the model from both space and time; (2) to understand temporal and spatial evolution characteristics of LSTs and explore the multiscale, abrupt and periodic characteristics of LSTs in the time series; and (3) to analyze the impact of LUCC on LSTs. The results of this study can provide a reference for hydrological simulation research in arid areas and drought monitoring and management for meteorological delays and other aspects.

Study Area
The Ebinur Lake Watershed (43 • 38 ~45 • 52 N and 79 • 53 ~85 • 02 E) is located in the arid area of Northwest China (Figure 1). The study area has a typical continental arid climate, characterized by dramatic annual temperature change, scarce precipitation, and strong potential evaporation [32,33], and the lake basin is surrounded by mountains on the south, west, and north sides, which block outside airflow [34]. Climate change easily affects the basin's water resources and ecosystems. The Ebinur Lake watershed is an important ecological barrier for environmental changes in the Junggar Basin in Xinjiang Uygur Autonomous Region (XUAR), which is a typical inland lake basin in an arid region [35]. Some studies found that the LUCC of the watershed has changed greatly [36]. However, the influence of LUCC on LST in the Ebinur Lake Watershed is not yet clear and needs to be further investigated. influence of LUCC on LST in the Ebinur Lake Watershed is not yet clear and needs to be further investigated.

Data Sources
The data used in this study can be described as follows: (1) The elevation data are based on the 90-m digital elevation model (DEM) of Shuttle Radar Topographic Mission (SRTM) data and were obtained from the Geospatial Data Cloud Platform of the Computer Network Information Center of the Chinese Academy of Sciences (http://www.gscloud.cn/accessed on 27 November 2021). (2) The global land cover data, with a resolution of 1 km, were obtained from the University of Maryland (http://app.earth-observer.org/data/basemaps/images/global/LandCover_512/Land-CoverUMD_512/LandCoverUMD_512.html/accessed on 27 November 2021), and the Terra and Aqua combined Moderate Resolution Imaging Spectroradiometer (MODIS) Land Cover Type (MCD12Q1) Version 6 data from the USGS EROS Center (http://doi.org/10.5067/MODIS/MCD12Q1.006/accessed on 27 November 2021), has a resolution of 500 m and uses two images (2001 and 2017). (3) The soil data were taken from the 0.5-degree soil data available for the VIC model released by the University of Washington (https://vic.readthedocs.io/en/master/Datasets/Datasets/accessed on 27 November 2021) and based on the World Soil Database (HWSD) 1:1 million Chinese soil dataset (v1.1). (4) The meteorological data were obtained from the daily dataset of surface climate data in China (V3.0), the time resolution is daily and the time scale is from 1 January 1960 to 31 December 2017. The data were obtained from the China Meteorological Data Network (http://data.cma.cn/accessed on 27 November 2021). (5) The verification data of LST were obtained from a combined Terra and Aqua MODIS LST and meteorological station data product for China from 2003 to 2017 [37] with a resolution of 5.6 km. The model resolution is 8 km*8 km, and the basin is divided into 767 grids ( Figure 2).

Data Sources
The data used in this study can be described as follows: (1)

VIC Model
The VIC hydrological model is a large-scale hydrological model with a variable infiltration capacity. The VIC model is an open-source large-scale distributed conceptual hy-

VIC Model
The VIC hydrological model is a large-scale hydrological model with a variable infiltration capacity. The VIC model is an open-source large-scale distributed conceptual hydrological model based on one kind of water and heat balance and physical dynamic mechanisms (VIC's homepage: http://vic.readthedocs.org/en/master/ accessed on 27 November 2021). Specific theoretical references include [38,39].
The energy balance calculation of the VIC-3L model mainly consists of four parts: LST, sensible heat flux, latent heat flux, and geothermal change. The other three parts are determined by LST, and the connecting factor of water balance and energy balance is latent heat flux. In an ideal situation, the energy balance calculation formula of each land cover type in the grid unit is as follows: where, R n is net radiation (W/m 2 ), H is sensible heat flux (W/m 2 ), ρ w L e E is latent heat flux (W/m 2 ), and G is surface heat flux (W/m 2 ). The net radiation is calculated by the following formula: where, α is the surface reflectance (%) of a certain land cover type, R S is the downward short-wave radiation (µm), ε is the specific emissivity(%), R L is the downward long-wave radiation (µm), and σ is the Stefan-Boltzmann constant (σ = 5.67 × 10 −8 W/(m 2 · • C 4 )). The calculation formula of sensible heat is as follows: where, Ts is the LST ( • C), Ta is the atmospheric temperature ( • C), and r h is the aerodynamic impedance (s/m). The surface heat flux is estimated by the heat of two layers of soil. For the upper layer (the first layer + the second layer), the total depth is assumed to be D 1 (m), then its calculation formula is: where, κ is the soil thermal conductivity(W/m· • C), T 1 is the temperature of the upper soil ( • C). Assuming that the depth of the underlying soil is D 2 (m) and the temperature at the bottom of the soil is constant, then: where, C S is the coefficient of soil heat conductivity (J/m 3 · • C), T + 1 and T − 1 is the soil temperature ( • C) at the beginning and the end of the period in the upper soil respectively, and T 2 is the constant temperature ( • C) of the lower soil. Currently, it is believed that κ and C S do not change with the change of soil moisture content. According to Equations (4) and (5), the calculation formula of surface heat flux can be deduced as follows: Combined with the above formula, the LST, sensible heat flux and surface heat flux can be calculated through iteration.

Morlet Continuous Wavelet
The long-term LST changes have obvious periodic characteristics, and there is less analysis of the multi-scale time characteristics between LST and meteorological elements. The wavelet analysis method has a powerful multi-scale resolution function, which can detect the response period of meteorological elements and LST, and can also show multiple change cycles hidden in the time series, reflecting the changing trend of the system within the same time range. Therefore, the wavelet analysis is more advantageous when studying the multi-scale characteristics of LST and meteorological elements. In this paper, we used the toolbox provided by Grinsted et al. [40], you can download it from this address: (http://www.glaciology.net/wavelet-coherence/ accessed on 27 November 2021).
The Morlet continuous wavelet function is a commonly used complex wavelet function. It has the advantages of rich phase information, good time aggregation, and highfrequency resolution. Therefore, it is widely used to study the correlation between two non-stationary sequences [41,42]. Therefore, this study selected the Morlet wavelet to analyze the degree of correlation between the LST and each meteorological element. Morlet wavelet is a complex wavelet, and its expressions in the time domain and frequency domain are: In the formula, ω 0 is a constant, and i is an imaginary number. Morlet wavelet is a single-frequency negative sine function under the Gaussian envelope. When using wavelets for feature extraction purposes the Morlet wavelet (with ω 0 = 6) is a good choice, since it provides a good balance between time and frequency localization.

Cross-Wavelet Transform (XWT)
Similar to the Fourier transform, according to the description of the transformed sequence energy in Perseval's theorem, that is, after the signal is transformed from the time domain to the frequency domain, the total energy remains unchanged, then the local time sequence x(t) The expression of Wavelet Power Spectrum (Wavelet Power Spectrum, WPS) is: The above-mentioned continuous wavelet transform is only a one-dimensional sequence of the wavelet transform, and this research mainly involves the correlation of twotime series. In the early research, it is necessary to perform a one-dimensional wavelet transform on each sequence separately, while Hudgins et al. [43], the XWT can not only directly analyze the correlation between two time series at different frequencies, but also reflect their phase information and local characteristics in the time domain and frequency domain. For any two series of time x(t) and y(t), the expression of the XWT is: In the formula, the wavelet transforms of the two time series x(t) and y(t) are W X t (ω) and W Y t (ω) respectively. The corresponding wavelet cross power (Cross Wavelet Power, XWP) is (XWP)xy = W XY t (ω) [44].

Wavelet Coherence (WTC)
WTC is to identify the possible relationship between two time series by finding the common points of two time series in different time scales and frequency ranges. It is a kind of enhanced time-frequency correlation. Recognition tool that can express the local coherence of the XWT in the time-frequency domain [45]. Even if they correspond to the low-energy region of the XWT, their correlation in the WTC spectrum may be significant. The expression of complex wavelet coherency (Complex Wavelet Coherency) ρ xy is: Remote Sens. 2021, 13, 4867 6 of 17 In the formula, S is a smooth spectrum operation in the time domain and the frequency domain. WTC is the absolute value of complex WTC. The mathematical expression is:

Accuracy Verification
The selection of the best model parameters is the main factor that affects the simulation accuracy of the VIC model. There are seven main parameters that need to be input to the VIC model: the storage volume curve index B; the maximum base flow velocity D m (m 3 /s); the nonlinear base flow velocity D s (m 3 /s); the soil moisture content W s (%); and the three layers of soil thickness d 1, d 2 and d 3 (m) during the nonlinear base flow. The selection of model parameters was based on traditional parameter calibration methods [46]. The above parameters were adjusted according to the empirical interval of the hydrological parameters to make the peak discharge value match the measured data. The above steps were repeated several times until a satisfactory result was obtained. The parameters are calibrated by repeatedly calculating the NSE and R 2 of the simulated and the measured values of runoff. The NSE and R 2 of the monthly runoff of Wenquan hydrological station during the verification period (2002-2017) are 0.44 and 0.47, respectively. The Final selected parameters are described as follows (Table 1):   (Figure 3d). It can be seen from the respective boxplots (Figure 3e,f) that the difference between the simulated and measured monthly data is not very large. The medians are 11.11 °C and 11.65 °C respectively. The measured values are slightly higher. The simulated and measured annual values have a large gap. They are 8.56 °C and 10.07 °C, the difference is 1.5 °C, our simulated values are underestimated.

Accuracy Verification in Terms of Space
In addition to verifying the LST over time, its spatial accuracy should be verified as well. Therefore, according to the remote sensing LST inversion dataset in China, the simulated LST based on the VIC model and its spatial distribution were compared and studied. However, through verification, it was found that the accuracy of the LST simulated by the VIC model for cold months is inaccurate, especially in January and December,

Accuracy Verification in Terms of Space
In addition to verifying the LST over time, its spatial accuracy should be verified as well. Therefore, according to the remote sensing LST inversion dataset in China, the simulated LST based on the VIC model and its spatial distribution were compared and studied. However, through verification, it was found that the accuracy of the LST simulated by the VIC model for cold months is inaccurate, especially in January and December, where the LST is even greater than 0 in high-altitude mountainous areas in winter. According to Koch's research [30], who found that the warmth of the VIC model in cold months is obvious, VIC has a distinct seasonality in its bias with a warm bias in winter (~3 • C) as opposed to a slight cool bias in summer (~−1 • C), which is related to an underestimation of the actual evapotranspiration and altitude. Our study area is located in an arid area, and the actual evapotranspiration in winter is small, which does not affect the accuracy of the LST simulated by the VIC model for areas other than high-altitude mountainous areas. To ensure the reasonableness of the data, the simulated LST values of the high-altitude mountainous area in November, January and December were deleted and converted to a no-data area. Taking 2017 as an example, the spatial resolutions were 0.05 • and 1/12 • respectively, as shown in Figure 4. It can be seen that the LST simulated by the VIC model has a high spatial consistency with the measured LST. The LST in the higher altitude area is significantly lower than that in other areas. As the time changes, the LST inside and outside the basin both "increase first and then decrease", and the simulated LST is slightly lower than the measured data. racy of the LST simulated by the VIC model for areas other than high-altitude mountain-ous areas. To ensure the reasonableness of the data, the simulated LST values of the highaltitude mountainous area in November, January and December were deleted and converted to a no-data area. Taking 2017 as an example, the spatial resolutions were 0.05° and 1/12° respectively, as shown in Figure 4. It can be seen that the LST simulated by the VIC model has a high spatial consistency with the measured LST. The LST in the higher altitude area is significantly lower than that in other areas. As the time changes, the LST inside and outside the basin both "increase first and then decrease", and the simulated LST is slightly lower than the measured data.
Of course, the accuracy needs to be improved, but for the data-scarce area, the simulation results can provide continuous and long-term serial data support, which has certain advantages.

Temporal and Spatial Variation Characteristics of the LST
The overall trend of the daily average LST in the study area was a unimodal distribution. The daily average LST ( Figure 5, left) began to increase at the end of February, reached the maximum in mid-July, began to decline at the end of August, and reached the minimum at the end of January. The peak in mid-July (July 13) was 25.92 °C, and the LST dropped rapidly to the lowest value in late January (January 22) at −18.99 °C. The daily average LST of the entire study area was 6.95 °C. Of course, the accuracy needs to be improved, but for the data-scarce area, the simulation results can provide continuous and long-term serial data support, which has certain advantages.

Temporal and Spatial Variation Characteristics of the LST
The overall trend of the daily average LST in the study area was a unimodal distribution. The daily average LST ( Figure 5, left) began to increase at the end of February, reached the maximum in mid-July, began to decline at the end of August, and reached the minimum at the end of January. The peak in mid-July (July 13) was 25.92 • C, and the LST dropped rapidly to the lowest value in late January (January 22) at −18.99 • C. The daily average LST of the entire study area was 6.95 • C.
The monthly average LST in the study area also showed a single-peak trend ( Figure 5  A histogram was used to represent the distribution of the annual average LST from 1960 to 2017 ( Figure 5, right). The figure shows that the average annual LST is on the rise, with a rising rate of approximately 0.0279 °C/a and an increase of approximately 1.62 °C in the past 58 years. The annual average LST is 6.40 °C. The annual average LST reached a trough in 1960, with a value of 4.54 °C; the figure clearly shows that the annual average LST after 2000 was significantly higher than that before 2000, which may be caused by increasing human activities, the industrial and agricultural activities in the Ebinur Lake Watershed. Based on the values of the daily LST simulated by the VIC model, the spatial distribution map of the inter-decadal and multi-year average LSTs from 1960 to 2017 was calculated. It can be seen that there are no obvious differences in the LST in the basin, but there are varying degrees of variation between different decadal periods.
The LST changes within the basin are not obvious, and the degree of change is different ( Figure 6). The main manifestation is the gradual decrease of the high LST areas in the east and the gradual increase of the central LST. The LST range changed from −11.48-9.04 °C in the 1960s to −11.14-9.56 °C in the 1980s, and then to −3.5-10.51 °C in 2010-2017, with the average values being 3.28 °C, 5.24 °C, and 6.79 °C, respectively. These results show that the LST in this area has an increasing trend. Although the range of high-value areas in the east has decreased, gradually shrinking to the low vegetation coverage area near Lake Ebinur, the LST value has increased, and the multi-year average LST in other regions has increased compared with the values in the 1960s. Therefore, the LST of the basin is increasing. Based on the values of the daily LST simulated by the VIC model, the spatial distribution map of the inter-decadal and multi-year average LSTs from 1960 to 2017 was calculated. It can be seen that there are no obvious differences in the LST in the basin, but there are varying degrees of variation between different decadal periods.
The LST changes within the basin are not obvious, and the degree of change is different ( Figure 6). The main manifestation is the gradual decrease of the high LST areas in the east and the gradual increase of the central LST. The LST range changed from −11.48-9.04 • C in the 1960s to −11.14-9.56 • C in the 1980s, and then to −3.5-10.51 • C in 2010-2017, with the average values being 3.28 • C, 5.24 • C, and 6.79 • C, respectively. These results show that the LST in this area has an increasing trend. Although the range of high-value areas in the east has decreased, gradually shrinking to the low vegetation coverage area near Lake Ebinur, the LST value has increased, and the multi-year average LST in other regions has increased compared with the values in the 1960s. Therefore, the LST of the basin is increasing.
Remote Sens. 2021, 13, x FOR PEER REVIEW 1 Figure 6. Interdecadal changes and annual mean LST spatial distribution in the Ebinur Lake shed. Figure 6. Interdecadal changes and annual mean LST spatial distribution in the Ebinur Lake Watershed.

Delayed Correlation Analysis of LST and Meteorological Elements
The XWT and WTC of the dual time series were applied to the coefficients after continuous wavelet transform, and the XWT spectrum between the LST of the Ebinur Lake Watershed and other meteorological elements was obtained ( Figure 7) and crossed with the WTC spectrum (Figure 8). The XWT and WTC can well reflect the "correlation" between changes in two different time series. The XWT generally reflects the strength of the "shared period" between sequences. The WTC generally reflects the consistency of the periodic "change trend" between sequences but does not directly reflect the intensity relationship of the change period. As shown in Figures 7 and 8 cycle of 1~1.5a (1964)(1965)(1966), and the two are negatively correlated; the other periodic feature is not obvious; the relative humidity and surface temperature have a significant area of 1~2.5a from 1964 to 1968 and a negative correlation, and the other two periodic features are not obvious; the sunshine hours and surface temperature have only one area where the periodic characteristics are not obvious. , and only the energy spectrum within this solid line needs to be considered. In order to eliminate the interference of the boundary effect, the thick solid line represents the key value that passes the significance test with 95% confidence; the arrow indicates the relative phase difference: → indicates that the two changes in phase are consistent, ← indicates that the two changes in phase are opposite, ↑ indicates that the phase change of the meteorological element is 90° ahead of the LST phase, and ↓ indicates that the phase change of the meteorological element is 90° behind the LST phase , and only the energy spectrum within this solid line needs to be considered. In order to eliminate the interference of the boundary effect, the thick solid line represents the key value that passes the significance test with 95% confidence; the arrow indicates the relative phase difference: → indicates that the two changes in phase are consistent, ← indicates that the two changes in phase are opposite, ↑ indicates that the phase change of the meteorological element is 90 • ahead of the LST phase, and ↓ indicates that the phase change of the meteorological element is 90 • behind the LST phase [47,48]). (The same below). [47,48]). (The same below).

The Impact of LUCC on LST
The statistics were classified according to the three-phase LUCC in this study. The specific changes in land cover types in the 1990s, 2001, and 2017 are shown in Figures 9  and 10, and the range distribution of the LST in the corresponding time period is shown in Figure 11.
The order of land cover types in the Ebinur Lake Watershed in the 1990s (Figure 9a) was open shrubland > grassland > bare ground > woodland grassland > crop land > water > woodland > evergreen needleleaf forest > mixed forest > urban and built-up > deciduous broadleaf forest. Among them (Figure 10), the largest area of open shrubland was 18,587.6 km 2 , accounting for 35.1% of the total area, followed by grassland, with an area of 18,476.3 km 2 , accounting for 34.9% of the total area. Bare ground accounted for 13.9% of the total area, and cropland accounted for 3.8% of the total area. Deciduous broadleaf forest occupied the smallest area, accounting for only 0.14% of the total area.  (1984-1990). The two showed a positive correlation between the average temperature and the LST; the annual precipitation and the LST had significant areas of 0.5~2a and 6~8a in 1964-1972 and 1984-1994, respectively. Both have an approximately negative correlation; relative humidity and annual precipitation are the same, and there are significant areas of 1~2a and 6~9a in 1964-1971 and 1978-1994, respectively, and the two also show a similar negative correlation. Sunshine hours and LST only existed in a significant area of 1~2a from 1964 to 1970, and the two are in an approximately positive correlation.
In the coherence spectrum of WTC, the average wind speed and surface temperature have a significant area greater than 15 years from 1980 to 1996, and the average wind speed in this area is ahead of the surface temperature, that is, 1/4a ahead of time; the average temperature and the surface temperature in the area with the 95% confidence interval accounted for more than 60%, indicating that the average temperature is the main controlling factor of the surface temperature, and the two existed in [1963][1964][1965][1966][1967][1968][1969][1970][1971][1972][1980][1981][1982][1983][1984][1985][1986][1987][1988][1989][1990][1991][1992][1993][1994][1995][1996] and 2004-2010 on the time scale of 1~4a. For significant periods, most areas of the two are approximately positively correlated. Among them, there may be abrupt changes around 1980 and 1996 on the time scale of 1~2a; that is, there may be a "high-low" transition of surface temperature, and there are still two significant cycles of~10a  and 5~7.5a (2000-2009); the annual precipitation and surface temperature only have a small cycle of 1~1.5a (1964)(1965)(1966), and the two are negatively correlated; the other periodic feature is not obvious; the relative humidity and surface temperature have a significant area of 1~2.5a from 1964 to 1968 and a negative correlation, and the other two periodic features are not obvious; the sunshine hours and surface temperature have only one area where the periodic characteristics are not obvious.

The Impact of LUCC on LST
The statistics were classified according to the three-phase LUCC in this study. The specific changes in land cover types in the 1990s, 2001, and 2017 are shown in Figures 9 and 10, and the range distribution of the LST in the corresponding time period is shown in Figure 11.
The order of land cover types in the Ebinur Lake Watershed in the 1990s (Figure 9a) was open shrubland > grassland > bare ground > woodland grassland > crop land > water > woodland > evergreen needleleaf forest > mixed forest > urban and built-up > deciduous broadleaf forest. Among them (Figure 10), the largest area of open shrubland was 18,587.6 km 2 , accounting for 35.1% of the total area, followed by grassland, with an area of 18,476.3 km 2 , accounting for 34.9% of the total area. Bare ground accounted for 13.9% of the total area, and cropland accounted for 3.8% of the total area. Deciduous broadleaf forest occupied the smallest area, accounting for only 0.14% of the total area. The order of land cover types in the Ebinur Lake Watershed in 2001 (Figure 9b) was grassland > non-vegetation land > cropland > water > woodland > urban and built-up > permanent wetland > evergreen needleleaf forest > deciduous broadleaf forest > mixed forest > woodland grassland > open shrubland > deciduous needleleaf forest. Among them (Figure 10), the area of grassland has increased, reaching the largest area of 34,714.9 km 2 and accounting for 65.6% of the total area, followed by non-vegetation land with an area of 11,904.9 km 2 , accounting for 22.5% of the total area. The two accounted for 88.1% of the total area. The area of cropland decreased to 4184.4 km 2 , accounting for only 7.9% of the total area. The area of open shrubland dropped sharply, most of them were transformed into grassland or non-vegetation land.
The order of land cover types in the Ebinur Lake Watershed in 2017 (Figure 9c) was grassland > non-vegetation land > cropland > water > woodland > urban and built-up > permanent wetland > evergreen needleleaf forest > deciduous broadleaf forest > mixed forest > woodland grassland > open shrubland > deciduous needleleaf forest. This is consistent with the ranking in 2001, and the area has not changed much. Grassland and nonvegetation land are also the main land cover types in the study area ( Figure 10).  In the 90s, the LST of grassland ranged from 5.2 to 8.4 °C, with a median of 6.3 °C; the LST of cropland ranged from 5.3 to 9.5 °C, with a median of 8.1 °C; the LST of nonvegetation land was in the range of 6.5-9.7 °C, with a median of 8.4 °C.
In 2001, the LST of grassland ranged from 5.7 to 8.5 °C, with a median of 6.6 °C; the The order of land cover types in the Ebinur Lake Watershed in 2001 (Figure 9b) was grassland > non-vegetation land > cropland > water > woodland > urban and built-up > permanent wetland > evergreen needleleaf forest > deciduous broadleaf forest > mixed forest > woodland grassland > open shrubland > deciduous needleleaf forest. Among them (Figure 10), the area of grassland has increased, reaching the largest area of 34,714.9 km 2 and accounting for 65.6% of the total area, followed by non-vegetation land with an area of 11,904.9 km 2 , accounting for 22.5% of the total area. The two accounted for 88.1% of the total area. The area of cropland decreased to 4184.4 km 2 , accounting for only 7.9% of the total area. The area of open shrubland dropped sharply, most of them were transformed into grassland or non-vegetation land.
The order of land cover types in the Ebinur Lake Watershed in 2017 (Figure 9c) was grassland > non-vegetation land > cropland > water > woodland > urban and built-up > permanent wetland > evergreen needleleaf forest > deciduous broadleaf forest > mixed forest > woodland grassland > open shrubland > deciduous needleleaf forest. This is consistent with the ranking in 2001, and the area has not changed much. Grassland and nonvegetation land are also the main land cover types in the study area ( Figure 10).  In the 90s, the LST of grassland ranged from 5.2 to 8.4 °C, with a median of 6.3 °C; the LST of cropland ranged from 5.3 to 9.5 °C, with a median of 8.1 °C; the LST of nonvegetation land was in the range of 6.5-9.7 °C, with a median of 8.4 °C.
In 2001, the LST of grassland ranged from 5.7 to 8.5 °C, with a median of 6.6 °C; the tation land was in the range of 6.7-9.6 °C, with a median of 7.8 °C.
During the period, the area of open shrubland, most of which was converted to grassland, cropland, and non-vegetation land, declined sharply; the area of grassland increased most obviously, and the LST was lower. Although the area of cropland was small, the LST was higher; compared with land types with higher vegetation coverage such as grassland, the LST of non-vegetation land was higher, which may have been because a large amount of solar radiation was absorbed by vegetation, lowering its LST.

Harmful Effects Caused by Changes in LST in Arid Regions
On the global scale, arid regions cover approximately one-third of the earth [49], and are one of the most sensitive regions to wind erosion phenomena [50]. In some recent research, a positive association was observed between the LST and wind erosion events [16,18]. The soil system has a synergistic relationship between various physical and chemical processes and environmental variables that affect its development [51]. Increasing the LST and decreasing the soil moisture, especially in the warm seasons over the arid regions, can be a factor in exacerbating soil erosion due to their effect on the salt accumulation on the soil surface [16,52]. Meanwhile, Different LUCC trends have different effects on the LST [53]. Heat sources (including bare land, urban and built-up land, and cropland) increase the LST, while cold sources (woodland, grassland, and water) reduce the LST [54]. Cropland and bare ground have a great impact on the regional LST, which is mainly reflected in the area and intensity of the impact [55]. This study found that between 1960 and 2017, the amount of open shrubs in the Ebinur Lake Watershed decreased sharply, and most of these areas were converted to grassland, cropland, and non-vegetation land. The vegetation coverage in the above-mentioned areas was low, and the LST was relatively high. This finding reasonably explains the variation in the LST of the Ebinur Lake Watershed from 1960 to 2017. However, the increased LST will lead to the reduction of regional soil moisture, vegetation degradation, and even increase the frequency of wind erosion and dust emission events. We believe this should attract the attention of local society.

Evaluation of VIC Model for LST Simulation
The LST is an important and complex hydrological state variable of land-atmosphere interactions. Land surface models (LSMs), such as the VIC model, are a key tool for enhancing the understanding of the process and providing predictions of the coupling between the terrestrial hydrosphere and its atmosphere. Distributed LSMs can be used to understand the hydrological status and flux of each grid cell, such as the LST. By comparing and analyzing the simulation accuracy of three different LSMs, Koch et al. [30] found The order of land cover types in the Ebinur Lake Watershed in 2001 (Figure 9b) was grassland > non-vegetation land > cropland > water > woodland > urban and built-up > permanent wetland > evergreen needleleaf forest > deciduous broadleaf forest > mixed forest > woodland grassland > open shrubland > deciduous needleleaf forest. Among them (Figure 10), the area of grassland has increased, reaching the largest area of 34,714.9 km 2 and accounting for 65.6% of the total area, followed by non-vegetation land with an area of 11,904.9 km 2 , accounting for 22.5% of the total area. The two accounted for 88.1% of the total area. The area of cropland decreased to 4184.4 km 2 , accounting for only 7.9% of the total area. The area of open shrubland dropped sharply, most of them were transformed into grassland or non-vegetation land.
The order of land cover types in the Ebinur Lake Watershed in 2017 (Figure 9c) was grassland > non-vegetation land > cropland > water > woodland > urban and built-up > permanent wetland > evergreen needleleaf forest > deciduous broadleaf forest > mixed forest > woodland grassland > open shrubland > deciduous needleleaf forest. This is consistent with the ranking in 2001, and the area has not changed much. Grassland and non-vegetation land are also the main land cover types in the study area ( Figure 10).
In the 90s, the LST of grassland ranged from 5.2 to 8.4 • C, with a median of 6.3 • C; the LST of cropland ranged from 5.3 to 9.5 • C, with a median of 8.1 • C; the LST of nonvegetation land was in the range of 6.5-9.7 • C, with a median of 8.4 • C.
In 2001, the LST of grassland ranged from 5.7 to 8.5 • C, with a median of 6.6 • C; the LST of cropland ranged from 7.9 to 9.7 • C, with a median of 8.7 • C; the LST of non-vegetation land was in the range of 6.6-9.7 • C, with a median of 8.1 • C.
In 2017, the LST of grassland ranged from 5.5 to 8.0 • C, with a median of 6.2 • C; the LST of cropland ranged from 7.1 to 8.7 • C, with a median of 8.3 • C; the LST of non-vegetation land was in the range of 6.7-9.6 • C, with a median of 7.8 • C.
During the period, the area of open shrubland, most of which was converted to grassland, cropland, and non-vegetation land, declined sharply; the area of grassland increased most obviously, and the LST was lower. Although the area of cropland was small, the LST was higher; compared with land types with higher vegetation coverage such as grassland, the LST of non-vegetation land was higher, which may have been because a large amount of solar radiation was absorbed by vegetation, lowering its LST.

Harmful Effects Caused by Changes in LST in Arid Regions
On the global scale, arid regions cover approximately one-third of the earth [49], and are one of the most sensitive regions to wind erosion phenomena [50]. In some recent research, a positive association was observed between the LST and wind erosion events [16,18]. The soil system has a synergistic relationship between various physical and chemical processes and environmental variables that affect its development [51]. Increasing the LST and decreasing the soil moisture, especially in the warm seasons over the arid regions, can be a factor in exacerbating soil erosion due to their effect on the salt accumulation on the soil surface [16,52]. Meanwhile, Different LUCC trends have different effects on the LST [53]. Heat sources (including bare land, urban and built-up land, and cropland) increase the LST, while cold sources (woodland, grassland, and water) reduce the LST [54]. Cropland and bare ground have a great impact on the regional LST, which is mainly reflected in the area and intensity of the impact [55]. This study found that between 1960 and 2017, the amount of open shrubs in the Ebinur Lake Watershed decreased sharply, and most of these areas were converted to grassland, cropland, and non-vegetation land. The vegetation coverage in the above-mentioned areas was low, and the LST was relatively high. This finding reasonably explains the variation in the LST of the Ebinur Lake Watershed from 1960 to 2017. However, the increased LST will lead to the reduction of regional soil moisture, vegetation degradation, and even increase the frequency of wind erosion and dust emission events. We believe this should attract the attention of local society.

Evaluation of VIC Model for LST Simulation
The LST is an important and complex hydrological state variable of land-atmosphere interactions. Land surface models (LSMs), such as the VIC model, are a key tool for enhancing the understanding of the process and providing predictions of the coupling between the terrestrial hydrosphere and its atmosphere. Distributed LSMs can be used to understand the hydrological status and flux of each grid cell, such as the LST. By comparing and analyzing the simulation accuracy of three different LSMs, Koch et al. [30] found that they all have certain spatial defects. The VIC has a distinct seasonality in its bias with a warm bias in winter (~3 • C) as opposed to a slight cool bias in summer (~−1 • C). The VIC model has a better simulation effect on the LST in warm months; in cold months, underestimation of the actual evapotranspiration will result in a high LST simulation. In this study, the same situation also appeared, but it may be because the study area is located in an arid area, that is, the actual evapotranspiration in autumn and winter, which are cold months, is inherently low, so except for high-altitude mountainous areas, the accuracy of the LST in cold months in plain areas is still good. However, the LST in high-altitude mountainous areas is obviously overestimated, only in November, December and January there is a big deviation in the mountainous area, and the LST in the areas close to the Tianshan Mountains is greater than zero in winter, which is unlikely to happen. This lack of spatial accuracy may be related to the occurrence of snow in mountainous areas, but further analysis is needed to study this in more detail.

Conclusions
This study used the VIC model to simulate the LST of the Ebinur Lake Watershed from 1960 to 2017, analyzed its temporal and spatial characteristics, conducted a multiscale crossanalysis of various meteorological elements and LST, analyzed the correlations between LUCC and LST, and obtained the following conclusions: 1.
Although the model simulates the LST well in time and the space verification results are relatively good, the LST simulation in the high-altitude area of the cold month is seriously overestimated, which may be related to the occurrence of snowfall or to the altitude. Further research is needed.

2.
The LST of the Ebinur Lake Watershed shows an overall increasing trend, and the annual average LST is higher in the central and eastern parts of the basin. On the temporal scale, the daily and monthly average LSTs showed unimodal trends. The interdecadal monthly changes are not obvious, and the monthly average LST from 2010 to 2017 fluctuates more than in other periods.

4.
From 1960 to 2017, the LUCC of the Ebinur Lake Watershed underwent major changes, and the reduction of open shrubs may have caused the LST increase in this area.