Changes in Lake Area in Response to Climatic Forcing in the Endorheic Hongjian Lake Basin, China

Endorheic lakes are key components of the water cycle and the ecological system in endorheic basins. The endorheic Hongjian Lake wetland is China’s national nature reserve for protecting the vulnerable species of Relict Gull. The Hongjian Lake, once China’s largest desert freshwater lake, has been suffering from severe shrinkage in the last two decades, yet the variations in the lake area and its responses to climate change are poorly understood due to a lack of in situ observations. In this study, using Landsat remote sensing images, the Modified Normalized Difference Water Index, and nonparametric tests, we obtained the Hongjian Lake area changes on the annual, seasonal, and quasi-monthly scales during 1988–2014, analyzed the corresponding variations of the six climatic factors in the Hongjian Lake Basin (HJLB) using satellite-based products, and investigated the multi-scale response characteristics of lake area to climatic forcing using correlation analysis. The results showed that the lake area decreased during 1988–2014, and this process can be divided into two sub-stages, namely the first slight increasing sub-phase in 1988–1999 and the second significant declining sub-phase in 2000–2014. The shifts in patterns of the seasonal cycle had three types: as the natural rhythm of the lake changes has been broken by intensive human activities since the late 1990s, the natural bimodal type I has obviously changed into non-natural bimodal type II and unimodal type III, featured by a declining peak in July–September. The climatic wet/dry regime on multi-scales during 1988–2014 in the HJLB was generally warming and drying, mainly reflected by the increase in temperature (T), arid index (AI) and evaporation (ET0, ETa), and the decrease in the precipitation (Pre) and actual water difference (AWD). There were large differences in the climatic factors at different time scales, especially in the wet and dry seasons. When the lagged effect, the cumulative effect, and the lagged and cumulative combined effect were gradually considered, the correlation coefficient significantly increased, and the direction of the correlation coefficient became coincident with common sense. The correlation analysis identified a lag period of approximately 1–3 years on an annual scale, and a lag period of approximately 1–3 months on a monthly scale. This study could provide a certain scientific reference for climate change detection, water resource management, and species habitat protection in the HJLB and similar endorheic basins or inland arid regions.


Introduction
Endorheic lakes or inland lakes refer to either permanent or seasonal lakes formed in an endorheic (hydrologically landlocked) basin, which is a closed drainage basin lacking any water outflow to

Study Area
The Hongjian Lake Basin (HJLB), with a longitude of 109 • 30 38" E−110 • 5 46" E and latitude of 38 • 51 10" N−39 • 20 20" N, is located in the Ordos endorheic basin in the Yellow River Basin, China [34,42]. The HJLB spans the Inner Mongolia Autonomous Region and Shaanxi Province, with four towns (i.e., Xinjie, Taigesumu, Zhongji, and Erlintu) within the HJLB (Figure 1). The Hongjian Lake, shaped like an irregular triangle, is located in the center of the HJLB. Seven seasonal rivers without in situ hydrometric stations flow into the central Hongjian Lake, namely the Maogaitu River, the Zhashake River, the Mudushili River, the Songdao River in the Inner Mongolia Autonomous Region, and the Qibosu River, the Erlintu River, and the Haolai River in Shaanxi Province ( Figure 1) [34,39,42]. Based on the hydrological analysis under the projected coordinate system projection of Beijing_1954_GK_Zone_19N, the area of HJLB is approximately 1513.10 km 2 , with an elevation from 1337 m to 1640 m [42]. The HJLB is a typical inland semi-arid and arid region, with an average annual precipitation of 393.0 mm in 1960-2014 and an average annual temperature from −10.0 • C in January to 7.4 • C in July [42]. The primary vegetation types in the HJLB are grassland, desert, woodland, cropland, and lake water surface [34,42,43].
The Hongjian Lake ( Figure 1) was officially recognized as a provincial scenic spot in 1995, built as a nature reserve at a county level in 1996, listed as an important wetland in the "Chinese Wetland Protection Motion Plan" in 2000, incorporated in the list of important wetlands in Shaanxi province and the list of important wetlands in China in 2008, upgraded to a national 4A-level scenic spot in 2012, became a provincial-level nature reserve in 2014 [34], and was officially recognized as a national-level nature reserve in 2018. The Hongjian Lake, as the locality with the largest colonies of Relict Gull since 2003, has been playing a key role in saving the Ordos Relict Gull's breeding population [32][33][34]. Relict Gull (Larus relictus) was listed as an endangered species in Annex A of the Convention on International Trade in Endangered Species of Wild Fauna and Flora (CITES) in 1997 [44], the Chinese red list of endangered animal species in 1998 [45], and as a "vulnerable" species on the Red List of the International Union for the Conservation of Nature (IUCN) since 2000 [32]. Relict Gull restrictedly breeds only on islands in saline and slightly saline lakes with fluctuating water-levels, and it is very sensitive to the changes in lake water level and lake area [32].

Remote Sensing Lake Area
In the poorly-gauged endorheic HJLB, satellite remote sensing images play a vital role in the rapid and accurate extraction of long-term lake areas. In this study, we used Landsat images to extract lake areas. First, we downloaded Landsat images with path/row numbers of 127/33 from the 1980s to 2014 from the United States Geological Survey's (USGS) remote sensing image database (http://earthexplorer.usgs.gov/). Because there were not many remote sensing images available in the early 1980s in HJLB, we chose to use level 1 Terrain-Corrected (L1T) Landsat TM image during 1988 −1999, Landsat TM/ETM+ images during 2000−2012, and Landsat OLI images in 2013 and 2014. In total, we chose about 324 remote sensing images, in which more than 300 images had clear skies or only slight cloud contamination around the lake regions. Of these images, 144 were used in 1988-1999 and 180 were used in 2000-2014. The detailed acquisition dates for the Landsat remote sensing images used for the HJLB are provided in the reference of Liang and Yan [34]. Second, based on the available Landsat TM/ETM+ and OLI images, we further applied the Modified Normalized

Remote Sensing Lake Area
In the poorly-gauged endorheic HJLB, satellite remote sensing images play a vital role in the rapid and accurate extraction of long-term lake areas. In this study, we used Landsat images to extract lake areas. First, we downloaded Landsat images with path/row numbers of 127/33 from the 1980s to 2014 from the United States Geological Survey's (USGS) remote sensing image database (http://earthexplorer.usgs.gov/). Because there were not many remote sensing images available in the early 1980s in HJLB, we chose to use level 1 Terrain-Corrected (L1T) Landsat TM image during 1988−1999, Landsat TM/ETM+ images during 2000−2012, and Landsat OLI images in 2013 and 2014. In total, we chose about 324 remote sensing images, in which more than 300 images had clear skies or only slight cloud contamination around the lake regions. Of these images, 144 were used in 1988-1999 and 180 were used in 2000-2014. The detailed acquisition dates for the Landsat remote sensing images used for the HJLB are provided in the reference of Liang and Yan [34]. Second, based on the available Landsat TM/ETM+ and OLI images, we further applied the Modified Normalized Difference Water Index (MNDWI) method proposed by Xu [46] to measure the water area of Hongjian Lake. The MNDWI method has been widely used to extract water bodies [34,[47][48][49], which is a modification of the Normalized Difference Water Index (NDWI) [50]. The MNDWI is expressed as where Green is a green band and MIR is a middle infrared band. Different MNDWI threshold values from 0.30 to 0.40 were tested, and an empirical threshold of 0.35 was finally selected to segment the MNDWI images into lake water area and non-water area in this study. The detailed three-step procedure used to measure the lake water area and the corresponding error analysis are provided in the reference of Liang and Yan [34].

Climatic Factors, Site's Meteorological Data, and Satellite-Based Products
Surface flow is scarce in endorheic lake basins which spatially concur with arid and semi-arid climate regions, with limited precipitation but high potential evaporation [3]. Precipitation is the most important source of water, and evaporation is the most important water consumption in endorheic basins. Temperature, as one of the most important factors for characterizing observed climate change [51], is the most important factor affecting evaporation. Therefore, precipitation, temperature, evaporation, and combined factors can be seen as important indicators of climate change in endorheic or inland basins. When considering the main sources of water resources in the closed HJLB, water consumption, major meteorological factors affecting water dissipation, as well as the availability of data with the same spatio-temporal resolution, we selected six indicators to characterize climate change in the HJLB. Basic information on the six climatic factors is shown in Table 1. In these indices, the precipitation (Pre) represents the source and supply of water resources, the temperature (T) is the key factor affecting snowmelt and evaporation, the potential evaporation (ET 0 ) represents evaporation capacity, the actual evaporation (ET a ) represents actual water consumption, the aridity index (AI) stands for characterization index of dryness, and the actual water difference (AWD) stands for the natural difference between water supply and consumption.
The analysis of changes in the mentioned above six climatic factors at monthly, seasonal, and annual scales requires sufficient meteorological data. In the poorly-gauged endorheic HJLB, when taking into account the spatial and temporal resolution of satellite-based products and their applicability at these scales, we chose the following two remote sensing products to obtain the precipitation (Pre) and actual evaporation (ET a ) in the HJLB ( Table 1). The Multi-Source Weighted-Ensemble Precipitation version 2.0 (MSWEP v2.0) [52] was determined to work well in capturing the monthly and annual scale precipitation characteristics in typical gauge-scarce basins in China [53,54]. The Global Land Evaporation Amsterdam Model version 3.2 (GLEAM v3.2) [55,56] achieved more consistent results with the monthly gauge observations in the Yellow River Basin and across China [57]. The MSWEP v2.0 and GLEAM v3.2 used in this study had the same spatial resolution of 0.25 • × 0.25 • . For the temperature (T) and potential evaporation (ET 0 ), we used the site's interpolation data. There are no national meteorological stations inside the basin, and only ten national meteorological stations are around the HJLB [42]. We obtained the daily meteorological data during 1980-2014 at the ten national stations from the China Meteorological Data Service Center (CMDC) [58] and calculated the potential evaporation (ET 0 ) at the ten stations using the Penman-Monteith equation recommended by the Food and Agriculture Organization (FAO) [59], interpolated the stations' ET 0 to grid data with a resolution of 0.25 • × 0.25 • using the thin plate spline method [60] in ANUSPLIN version 4.4 software (Australian National University, Canberra, Australia), and then extracted the targeted ET 0 by the mask of HJLB. Similarly, the average temperature (T) in the HJLB is also obtained by the same interpolating method on the basis of the ten stations' data.
It is worth noting that, in order to better analyze the lagged response characteristics of the lake area to climatic factors, the time range of the meteorological factors  is longer than that of lake data . Furthermore, in this study, for the seasonal time scale, according to the actual conditions of rainfall and temperature in the HJLB, as well as the continuity of hydrologic time, the seasons are divided into the wet season (May−October) and the dry season (November−April).

Basic Statistical Measures
In this study, to describe the change characteristics of the lake area and the related climatic factors from the four aspects of mean, volatility, trend, and significance, we selected the mean value (X), standard deviation (SD), coefficient of variation (CV), Sen's slope of the nonparametric linear trend (S), and Mann-Kendall's Z as basic statistics. In particular, considering the different units among the lake area and the climatic factors, we used the CV that has the advantage of comparing distributions obtained with different units to better compare the variability [61,62]. We used the IBM SPSS Statistics Version 19 software to calculate X, SD, and CV. According to the definitions of S and Z and their calculation formulas (Section 2.4.2), we used the MATLAB R2017b (9.3.0.713579) software (The MathWorks, Inc., Natick, MA, USA) to calculate the S and Z.

Trend Analysis
We applied the nonparametric Sen's slope estimator method [63] to obtain the slope of linear trends in the lake area and the climatic factors at annual, seasonal, and monthly scales during 1988-2014 and the two sub-periods (1988-1999 and 2000-2014). Sen's slope estimator has been widely used in hydro-meteorological time series [64,65] due to its robustness for non-normally distributed data, and it is far less influenced by outlier and missing data. For the given data set X i (i = 1, 2, ..., n), the slope is first computed by where X j and X k are data values at times j and k (j > k), respectively. The n is the data set length, which is different at different time scales in different periods. For example, in the study period of 1988-2014, n equals 27 at the annual and seasonal scale, and 324 at the monthly scale, respectively. The Sen's slope (S), defined as the median of these n values of S i , is finally estimated as Then, we used the nonparametric Mann-Kendall test [66,67] to detect the significance of trends. The Mann-Kendall test is a widely used method for trend detection due to its robustness for non-normally distributed data [68][69][70][71]. For the same data set X i (i = 1, 2, ..., n), their corresponding ranks are R 1 , R 2 . . . R n , and the Mann-Kendall test statistic (Rs) is denoted as [64] where the symbolic function is defined as Mann and Kendall [66,67] have shown that when the number of data set n ≥ 8, the variance of the Mann-Kendall test statistic is calculated by where m is the number of tied groups and t i is the numbers of ties of extent i. The standardized test statistic Z is calculated as A positive Z indicates increasing trends, while a negative Z indicates decreasing trends in the time series. When |Z| > Z 1−α/2 , the null hypothesis is rejected and a significant trend exists in the time series. Z 1−α/2 is the critical value of Z from the standard normal distribution table. In this study, we chose α = 0.05 and α = 0.01 as the significance levels. The critical value of Z 1−α/2 is 1.960 at the 5% significance level and 2.576 at the 1% significance level.

Abrupt Change Point Detection
We applied the non-parametric Pettitt test [72] to detect the turning points in the series of the lake area and climatic factors. For the same given time series of X 1 , X 2 , . . . , X n , the Pettitt test sets the series as two samples X 1 , X 2 . . . , X t and X t+1 , X t+2 , . . . , X n , and the corresponding Mann-Whitney statistic U t,n is expressed as where the symbolic function sign(X t −X j ) is defined as where n is the data set length. Then, the statistic K t,n and the significance level are defined as When p < 0.05, a significant turning point exists at the 5% significance level. According to Equations (8)-(11), we used the MATLAB R2017b (9.3.0.713579) software to calculate the corresponding results.

Cumulative Anomaly Analysis
In this study, to analyze the cumulative effects of climatic factors on lake area changes, we used the cumulative anomaly analysis (CAA) to obtain the cumulative anomaly series of the lake area and the climatic factors on annual and monthly scales. The cumulative anomaly (CA) is calculated using the following formula [73,74]: where X i is the time series and X is the mean value of all the time-series data.

Spearman Rank Correlation Analysis and Design
We used the Spearman rank correlation analysis to measure the strength of the relationships between the lake area and the climatic factors. The Spearman rank correlation analysis [75] is essentially a nonparametric correlation test to detect nonlinear (monotonic) trends between two variables, with the strengths that it does not require a particular data distribution and is not sensitive to outliers. Spearman's rank correlation coefficient rho (ρ) is a nonparametric correlation coefficient ranging between −1.0 (perfect negative relationship) and 1.0 (perfect positive relationship). For a set of n sample pairs (X i , Y i ), the data for the two variables are ranked independently among themselves. We denote the ranks of X i as R X i , and ranks of Y i as R Y i , and ties in X i or Y i are assigned average ranks. The Spearman's rank correlation coefficient rho (ρ) is calculated as [75][76][77][78] We applied Ramsey's table of critical values for the Spearman's rank correlation coefficient rho (ρ) to test significance [76]. Ramsey's table [76] gives the exact critical values of ρ for 3 ≤ n ≤ 18 and Edgeworth approximations for 19 ≤ n ≤ 100. For example, in this study, for n = 27 (i.e., the sample pairs at an annual scale during 1988-2014), the critical value is 0.383 for a significance level of α = 0.05 and 0.492 for α = 0.01 (two-tailed); when n = 12 during 1988-1999, the critical value is 0.587 for α = 0.05 and 0.727 for α = 0.01 (two-tailed); when n = 15 during 2000-2014, the critical value is 0.521 for α = 0.05 and 0.654 for α = 0.01 (two-tailed). For n > 100, according to Nijsse and Zar [79,80], the following t-test can be used to test the significance of rho (ρ): where Z sp is the Student's t-distribution with (n − 2) degree of freedom. For the significance level of α = 0.05 and α = 0.01 (two-tailed), the corresponding critical value of the Student's t-distribution is defined as t (n−2, α) . In this study, for n = 324 (i.e., the sample pairs at a monthly scale during 1988-2014), the critical value is 1.967 at α = 0.05 and 2.591 at α = 0.01 (two-tailed); when n = 144 during 1988-1999, the critical value is 1.977 at α = 0.05 and 2.610 at α = 0.01 (two-tailed); when n = 180 during 2000-2014, the critical value is 1.973 at α = 0.05 and 2.604 at α = 0.01 (two-tailed). If Z sp > t (n−2, α) , the Spearman's rank correlation coefficient rho (ρ) is significant. We used the MATLAB R2017b (9.3.0.713579) software to calculate the ρ and obtain its significance (two-tailed), and validated the results using the IBM SPSS Statistics Version 19 software.
In order to better gradually investigate the lagged response characteristics of the lake area to climate factors and the cumulative effect of climate factors on lake area changes, we designed a series of progressive correlation analysis schemes for the data pairs (lake area vs. climatic factor), including four cases (Table 2). Case 1 is the original correlation analysis of the original data pairs, which is completely consistent in time (e.g., area in January vs. precipitation in January); case 2 is the lagged correlation analysis (e.g., the area in February (one month lag) vs. precipitation in January); case 3 is the cumulative correlation analysis (e.g., the cumulative anomaly of area in January-December vs. the cumulative anomaly of precipitation in January-December); case 4 is the lagged and cumulative correlation analysis, which is used first to perform the lag processing on the original area data (e.g., the area in February (one month lag) vs. precipitation in January), and then to carry out the cumulative anomaly processing on the lagged area sequence and the original precipitation sequence, and finally to do a correlation analysis of the data pairs of the transformed cumulative series. Table 2. Schemes of the Spearman correlation analysis designed in the Hongjian Lake Basin.

Cases
Name Time Series (Lake Area vs. Climatic Factor)

Case 1 Original correlation analysis
Original lake area series vs. original climatic factor series Case 2 Lagged correlation analysis Lagged lake area series vs. original climatic factor series Case 3 Cumulative correlation analysis Cumulative anomaly series of lake area vs. cumulative anomaly series of climatic factor

Case 4 Lagged and cumulative correlation analysis
Cumulative anomaly series of lagged lake area vs. cumulative anomaly series of climatic factor

Variations in Lake Area at Multiple Time Scales
We sequentially analyzed the changes in the average lake area at annual, seasonal, monthly, and intra-annual scales, plotted the process line for the change in the lake area ( Figure 2), and calculated the basic statistics (Table 3) at each scale to analyze its variations and characteristics.  Second, Figure 2b and Table 3 show the variations in the lake area in the wet season (May-October) and the dry season (November to April in the next year). During 1988-2014, whether in the wet or dry season, the overall change process and significant declining trend of the lake area were similar to that at the annual scale. The total significant downward trends were −0.99 km 2 /a in the wet season and −0.97 km 2 /a in the dry season. However, in the two sub-periods, there were still slight   First, Figure 2a and Table 3 show the variations in the annual lake area and its basic statistics. During 1988−2014, the Hongjian Lake presented a large fluctuation. The mean value (X), standard deviation (SD), and coefficient of variation (CV) of the annual lake area were 44.86 km 2 , 7.86 km 2 , and 0.18, respectively. The annual lake area had a significant declining trend (S) of −1.01 km 2 /year. When the Mann-Whitney statistic U t,n in the Pettitt test was examined (blue solid line in Figure 2), the maximum U t,n was around the year 2000, and its calculated p-value was less than 0.05 (Equation (11)), which meant that there was one significant abrupt change point around the year 2000 examined by the Pettitt test. Therefore, the whole change process in the lake area could be divided into two distinct sub-phases, namely the slight ascending sub-phase with a small fluctuation in 1988-1999 and the significant monotonic declining sub-phase in 2000-2014. Specifically, in the first sub-period of 1988-1999, the X, SD, CV, S, and Z of the annual lake area were 52.59 km 2 , 1.79 km 2 , 0.03, 0.07 km 2 /year, and 0.34 (not significant), respectively. In the second sub-period of 2000-2014, the average annual lake area (X) dramatically decreased to 38.68 km 2 , with a reduction ratio of 26.5% compared with that in 1988−1999. The volatility and change rate of the annual lake area during 2000-2014 increased, and the corresponding SD, CV, and S of the annual lake area changed to 4.48 km 2 , 0.12, and −1.00 km 2 /year, respectively. The Z was equal to -5.05, which was significant at the 1% significance level.
Second, Figure 2b and Table 3 show the variations in the lake area in the wet season (May-October) and the dry season (November to April in the next year). During 1988-2014, whether in the wet or dry season, the overall change process and significant declining trend of the lake area were similar to that at the annual scale. The total significant downward trends were −0.99 km 2 /a in the wet season and −0.97 km 2 /a in the dry season. However, in the two sub-periods, there were still slight differences. In the first sub-period (1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999), the average lake area in the wet season (52.75 km 2 ) was a little larger than that in the dry season (52.58 km 2 ), but the ascending trend in the wet season (0.03 km 2 /a) was smaller than that in the dry season (0.10 km 2 /a). In the second sub-period (2000−2014), the X, SD, CV, and significant S of the lake area in the wet season were 38.65 km 2 , 4.42 km 2 , 0.11, and −0.97 km 2 /year, respectively; the X, SD, CV, and significant S in the dry season changed to 39.06 km 2 , 4.58 km 2 , 0.12, and −1.00 km 2 /year, respectively. The absolute values of all the basic statistics in the wet season were smaller than that in the dry season.
Third, Figure 2c and Table 3 show the variations in the lake area at a quasi-monthly scale. The lake area on the quasi-monthly scale had more subtle changes, which appeared to be more substantial, informative and winding on the variation process line. During 1988-2014, the X, SD, CV, and significant trend (S) of the lake area on the quasi-monthly scale were 44.86 km 2 , 7.78 km 2 , 0.17, and −0.082 km 2 /month (i.e., −0.99 km 2 /year), respectively. In the first sub-period of 1988-1999, the lake area on the quasi-monthly scale belonged to a high-value and low-oscillation phase, namely with a high X (52.59 km 2 ) and low CV (0.04); while the second sub-period of 2000−2014 became a low-value and high-oscillation phase, namely with a low X (38.68 km 2 ) and high CV (0.11).
Lastly, the lake area at the intra-annual scale is shown in Figure 3 and Table 4. The intra-annual variations in the lake area during 1988−2014 had two peaks because there were two types of natural recharge in the lake basin, namely snowmelt runoff recharge in March−April and rainfall recharge in July−October. In the sub-period of 1988-1999, the intra-annual variations in the lake area also had two peaks, while in the second sub-period of 2000-2014, there was only one peak left in March-April (Figure 3a). When we took a closer look at the intra-annual changes each year, we found that the intra-annual variations in the lake area could be further divided into three types, namely type I, II, and III. Type I was a bimodal type, having one low peak on the left in March or April and one high peak on the right in July, August, or September. Type II was a bimodal type, with one high peak on the left in March or April and one low peak on the right in August, September, or October. Type III belonged to the unimodal type, with one single peak on the left in March or April and declining on the right. The years of type I were mainly concentrated before 2000, as well as the years of heavy precipitation in 2002 and 2003. The years of type II and III were concentrated after 2000, as well as some years of less precipitation (i.e., 1989, 1990, 1993, 1997, and 1999) (Table 4). For example, the intra-annual variations in the lake area in 1988 and 1991 belonged to type I, and the process of change exhibited a remarkable bimodal type: there was one low peak in March and one high peak in September in 1988, and there was one low peak in April and one high peak in July in 1991 (Figure 3b). The intra-annual variations in the lake area in 2004 and 2005 belonged to type II: there was one high peak in March or April and one low peak in August-October (Figure 3c). The intra-annual variations in the lake area in 1988 and 1991 belonged to type III: there was one peak in April, and then there was a stepped downward trend with no obvious peak after April (Figure 3d). statistics in the wet season were smaller than that in the dry season.
Third, Figure 2c and Table 3 show the variations in the lake area at a quasi-monthly scale. The lake area on the quasi-monthly scale had more subtle changes, which appeared to be more substantial, informative and winding on the variation process line. During 1988-2014, the , SD, CV, and significant trend (S) of the lake area on the quasi-monthly scale were 44.86 km 2 , 7.78 km 2 , 0.17, and − 0.082 km 2 /month (i.e. −0.99 km 2 /year), respectively. In the first sub-period of 1988-1999, the lake area on the quasi-monthly scale belonged to a high-value and low-oscillation phase, namely with a high (52.59 km 2 ) and low CV (0.04); while the second sub-period of 2000−2014 became a low-value and high-oscillation phase, namely with a low (38.68 km 2 ) and high CV (0.11).

Figure 3. (a)
The intra-annual variations in the Hongjian Lake area during different periods; (b) the type Ⅰ (bimodal: low peak on the left, high peak on the right), (c) type Ⅱ (bimodal: high peak on the left, low peak on the right), and (d) type Ⅲ (unimodal: single peak on the left, declining on the right) of intra-annual variations in some representative years.
Lastly, the lake area at the intra-annual scale is shown in Figure 3 and Table 4. The intra-annual variations in the lake area during 1988−2014 had two peaks because there were two types of natural recharge in the lake basin, namely snowmelt runoff recharge in March−April and rainfall recharge in July−October. In the sub-period of 1988-1999, the intra-annual variations in the lake area also had two peaks, while in the second sub-period of 2000-2014, there was only one peak left in March-April (Figure 3a). When we took a closer look at the intra-annual changes each year, we found that the intra-annual variations in the lake area could be further divided into three types, namely type Ⅰ, Ⅱ, and Ⅲ . Type Ⅰ was a bimodal type, having one low peak on the left in March or April and one high peak on the right in July, August, or September. Type Ⅱ was a bimodal type, with one high peak on the left in March or April and one low peak on the right in August, September, or October. Type Ⅲ belonged to the unimodal type, with one single peak on the left in March or April and Figure 3. (a) The intra-annual variations in the Hongjian Lake area during different periods; (b) the type I (bimodal: low peak on the left, high peak on the right), (c) type II (bimodal: high peak on the left, low peak on the right), and (d) type III (unimodal: single peak on the left, declining on the right) of intra-annual variations in some representative years. Table 4. Types of intra-annual variations in the Hongjian Lake area during 1988-2014.

Types a Characteristics Years
Type I Bimodal: low peak on the left, high peak on the right 1988,1991,1992,1994,1995,1996,1998,2002,2003 Type II Bimodal: high peak on the left, low peak on the right

Variations in Climatic Factors at Multiple Time Scales
Similarly, we sequentially analyzed the corresponding variation characteristics of the six climatic factors closely related to lake area changes at annual, seasonal, monthly, and intra-annual scales. The turning points of these factors' variation were different from that of the lake area. The turning point of precipitation (Pre) and actual evaporation (ET a ) was in the year 2001, and the turning point of temperature (T), potential evaporation (ET 0 ), aridity index (AI) and actual water difference (AWD) was in the year 1996 ( Figure 4). However, to be consistent with the analysis of changes in the lake area, we had to use the same sub-periods of the lake area change as the dividing line of the climatic factor change.

Correlation between Lake Area and Climatic Factors
The correlation analysis results between the lake area and the six climatic factors at different scales (i.e., annual, wet season, dry season, and monthly) during 1988−2014 and the two sub-periods (1988−1999 and 2000−2014) are shown in Table 5. The Spearman correlation coefficients at different scales showed large differences. Regarding the magnitude of the correlation coefficient, the correlation coefficients at all scales were very low. The overall range of the correlation coefficients was from −0.50 to 0.43, and the absolute value of most correlation coefficients was less than 0.30. Even the absolute value of most correlation coefficients on the monthly scale was less than 0.10; no matter in which period or sub-period, the absolute values of the correlation coefficients at the annual scale and the wet season scale were generally higher than at the dry season scale and the monthly scale. Regarding the directionality of the correlation coefficient, regardless of in the entire study period or its two sub-periods, there was always a certain factor in the direction of the correlation coefficient that contradicts common sense (italic numbers labeled in Table 5), showing the opposite direction to common sense. For example, as the most important water resources for the lake, the precipitation (Pre) should be positively related to the lake, so the negative relationship is abnormal. This reversal of common sense was particularly evident in the second sub-period (2000−2014). When further considering the lagged response, the lagged correlation analysis results between the lagged lake area and the six climatic factors at annual and monthly scales are shown in Table 6. In general, the lagged correlation coefficient was slightly increased compared to the normal (no-lagged) correlation coefficient, especially considering the lag period of nearly 3 years and nearly 3 months. Moreover, the directionality of some correlation coefficients that violate common sense had also changed. For example, the correlation coefficient between the lake area and Pre at the annual scale during 1988−2014 was −0.05, and the 3-year-lagged correlation coefficient was 0.17. It is worth noting that the lagged correlation coefficients at the monthly scale were still very low; in the second sub-period (2000−2014), the directionality of the correlation coefficient for most factors still appeared to be contrary to common sense.
When further considering the cumulative effect, the cumulative correlation analysis results between the cumulative lake area and the six cumulative climatic factors at multi scales are shown in Table 7. In general, the magnitude and significance of the cumulative correlation coefficient were increased compared to the normal (non-cumulative) correlation coefficient, especially for T, ET 0 , ET a , and AWD. For example, in the first sub-period of 1988−1999, the cumulative correlation coefficients for Pre, T, ET 0 , ET a , AI, and AWD at an annual scale were 0.30, −0.69 (significant), −0.52, −0.99 (significant), −0.45, and 0.80 (significant), respectively. Moreover, the cumulative correlation coefficients of all indicators on the monthly scale increased obviously. For example, during 1988−1999, the significant cumulative correlation coefficients for Pre, T, ET 0 , ET a , AI, and AWD at the monthly scale were 0.43, −0.53, −0.31, −0.81, −0.19, and 0.87, respectively. It is worth noting that the directionality of cumulative correlation coefficients for Pre and AI at all scales during 1988−2014 and on the annual, wet season, and monthly scale during 2000−2014 still appeared to be contrary to common sense. Table 6. Spearman correlation coefficient for analysis to the lagged response of lake area and climatic factors on annual and monthly scales. When further jointly considering the lagged effect and the cumulative effect, the integrated lagged-cumulative correlation analysis results between the integrated lake area and the six cumulative climatic factors at annual and monthly scales are shown in Table 8. The lagged-cumulative correlation coefficient at annual and monthly scales further improved in size, directionality, and significance. Taking the first sub-period with the most significant improvement effect as an example, during 1988−1999, the significant 2-year-lagged-cumulative correlation coefficients for Pre, T, ET 0 , ET a , AI, and AWD at the annual scale were 0.94, −0.98, −0.82, −0.70, −0.94, and 0.94, respectively. The significant 3-month-lagged-cumulative correlation coefficients for Pre, T, ET 0 , ET a , AI, and AWD at a monthly scale were 0.49, −0.60, −0.39, −0.79, −0.31, and 0.87, respectively. It is still worth noting that the directionality of lagged-cumulative correlation coefficients for Pre and AI at a monthly scale during 1988−2014 and on a monthly scale during 2000−2014 still slightly appeared to be contrary to common sense. Table 7. Spearman correlation coefficients of cumulative anomaly sequence between lake area and climatic factors at different scales during 1988−2014 and the two sub-periods a .    Overall, according to the corresponding results from the Tables 5-8, when we gradually considered the four cases or scenarios of progressive Spearman correlation analysis schemes (see in Table 2), namely case 1 (original correlation analysis), case 2 (lagged correlation analysis), case 3 (cumulative correlation analysis), and case 4 (lagged and cumulative correlation analysis), the lagged effect and cumulative effect were gradually reflected and displayed by the increase of the correlation coefficient, the reasonable change in the direction of the correlation coefficient, and the significant improvement of the correlation coefficient.

Lake Shrinkage in This Study and Other Endorheic Basins
Due to the influence of current climate warming and intensified human activities, whether in the Ordos endorheic (inflow) region where our study area is located or the larger inland area in northwest China and other famous inland regions worldwide, severe lake shrinkage has been encountered in recent decades. For example, in the Ordos endorheic region, Yan et al. [47] studied the changes of the Bojiang Lake area using the Landsat remote sensing data and revealed that the average lake area significantly declined from 8.86 km 2 in 1980-1998 to 3.50 km 2 in 1999-2013, with average change rates of -60.5% and a changing trend of −0.75 km 2 /year in 1999−2013. Liang and Yan [34] studied the changes in the annual Hongjian Lake area in the Ordos endorheic region using Landsat images and revealed that the average change rates and trends during 1999-2015, compared with those in 1988-1998, were -26.5% and -1.04 km 2 /year, respectively. In northwest China, Bai et al. [13] applied Landsat imageries to find that Bosten Lake, as the largest inland freshwater lake in China, declined from 1,055.83 km 2 in 1975 to 958.90 km 2 in 2007, with an average change rate of −9.18% and a changing trend of −2.94 km 2 /year. Guo et al. [19] also revealed that the area of Bosten Lake decreased by −348 km 2 from 2003 to 2010, with an average declining trend of −43.50 km 2 /year. In Central Asia, Bai et al. [13] used Landsat imageries to reveal the lake area of the Aral Sea, once the world's fourth-largest lake, which decreased from 59,262.35 km 2 in 1975 to 14,400.15 km 2 in 2007, with a change rate of −75.70% and average change trend of −1359.46 km 2 /year. Chaudhari et al. [15] used Landsat imageries to reveal that Urmia Lake in northwestern Iran, the largest lake in the Middle East and once the second-largest saline lake worldwide, decreased from approximately 5000 km 2 in 1987 to 1150 km 2 in 2016, with a shrinkage rate of −77% and average declining trend of −128.33 km 2 /year. In our study, Hongjian Lake declined on average by 52.99 km 2 annually in 1988-1999 to 38.68 km 2 in 2000-2014, with a declining ratio of -26.5%, and by −1.00 km 2 /year in 2000-2014. Therefore, the shrinkage rate and trend in our study were within the scope of these cited research results.
Comparing these studies, we find that lake wetland shrinkage is a relatively common phenomenon, especially in endorheic basins which spatially concur with arid and semi-arid climate regions. These endorheic lakes have suffered severe shrinkage, but there are still significant differences in the size and rate of change in the lake area due to differences in remote sensing data sources, data time ranges, trend estimation method, regional climatic conditions, and the anthropogenic disturbance intensity. These differences reflect the complexity and regional differences of inland lake changes. These regional differences and their hidden mechanisms are worthy of in-depth comparative research in the future. The issues regarding the lagged response and cumulative effect initially revealed in this study are worthy of further exploration. A comprehensive comparison of the lagged response characteristics and cumulative effects between different endorheic basins or regions worldwide will be very interesting.
Whether in the research cited in this study or in other related research that is not cited in this study, the nonparametric statistical approach is the preferred methodology for studies on the variations of the lake itself, climatic and anthropogenic factors, and the relationships between the lake and the factors. Nonparametric statistical methods have the advantages that they require few assumptions about the underlying populations, are often easier to apply than the normal theory counterparts, enable the user to obtain exact p-values for tests, have robustness for non-normally distributed data, and are insensitive to outlier and missing data [78]. The combined use of non-parametric methods has a certain reference significance for similar research and is highly recommended, such as the Sen's slope estimator method [63], the Mann-Kendall test [66,67], the Pettitt test [72], and the Spearman rank correlation test [75] used in this study. In addition, the design ideas and design schemes of the correlation analysis in this study (Section 2.4.5 and Table 2) play a key role in investigating the lagged response effect and cumulative effect, which has a certain reference significance for other related research on lake wetlands.

Implications for Lake Wetland Protection in Endorheic Basins
First, we reveal that both for the continuous linear decline of the annual lake area or the morphological change in the curve of intra-annual variation process, it is fully explained that the natural rhythm of the lake change has been broken, and the actual changes of lake area are not well explained by climatic forcing, which indicates that the lake area change is impacted more by the intensive human activities. The annual lake area has been decreasing year by year since 1999 or 2000, and this decline has almost continued to decline linearly, which is totally incompatible with the natural rhythm of climatic dry-wet conditions. Before 2000, the curve of the annual variation of the lake area clearly showed a "double peak" pattern or bimodal type, which is basically consistent with the natural intra-annual rhythm of climate change; however, this natural "double peak" pattern was obviously artificially broken after 2000. The "double peak" form has evolved into a "single peak" form or unimodal type, and even shows a linear decline pattern that decreases month by month. This pattern change will only recover to a certain degree when the precipitation is much larger than the average. For example, in 2002, 2003, 2012, and 2013, the precipitation in these years was  mm larger than the mean value. Liang et al. [34] showed sharp increasing variations in the annual human activities during 1988-2015, using the seven representative indicators (i.e., area of land use types, Normalized Differential Vegetation Index, human population, sheep population, gross industrial output value, hydraulic engineering construction, and tourism development), and revealed that the human activities had a more important influence on the lake area change after 1999. These studies have jointly demonstrated that, since the late 1990s or early 21st century, human activities have become the dominant factor in the change of the lake area. However, due to the lack of data on human activities on a monthly scale, we are unable to determine the response characteristics of the lake area to human activities. It is urgent to collect more detailed human activity data from multiple sources, and we recommend that the monitoring of water use data for human activities must be strengthened to provide available data support for subsequent research. Under the influence of such strong human activities, reasonably controlling human activities to reduce human water use in the basin has become a new core task of water resource management.
Second, to make matters worse, climate change has also cast a shadow over the protection of the endorheic lake wetland. Our results showed that, during 1988−2014, whether on the annual scale, dry and wet season scale, or monthly scale, the average climate of the basin was generally warming and drying, mainly reflected by the increase in temperature (T), arid index (AI) and evaporation (ET 0 , ET a ), and the decrease in the precipitation (Pre) and actual water difference (AWD) (Figure 4). When we checked this change on the intra-annual scale (2000−2014 vs. 1988−1999), we further found that Pre and AWD mainly reduced in March−May and July−August, which coincides with the peak of the snowmelt recharge in March and rainfall recharge in July−August of the lake, while T and ET a increased in most months except December. This climatic warming and drying have led to water deficits in the lake basin, and the lake area has accordingly shrunk on the annual scale, dry and wet season scale, and monthly scale since the late 1990s or early 21st century. We also found that lake changes have a lagging response to climate change, with a lag period of about 1−3 years on the annual scale and 1−3 months on a monthly scale. This means that the size of the lake area in a certain year is significantly affected by the climatic dry−wet regime of the previous three years, and the size of the lake area in a certain month is obviously affected by the climatic dry−wet conditions for the past three months. Hence, for a natural closed endorheic lake basin like the Hongjian Lake, without the supply of other external water sources, the change in the lake area is largely controlled by changes in the climatic conditions within the basin itself-especially the climatic dry−wet regime in recent months and in recent years. This high susceptibility to climate change is an inherent property of closed basins, and relevant authorities must clearly recognize the fundamental water scarcity and the natural vulnerability of water resources.

Conclusions
In this study, using Landsat remote sensing images and the Modified Normalized Difference Water Index, we obtained the Hongjian Lake area changes on annual, seasonal, and quasi-monthly scales during 1988−2014, analyzed the corresponding variations of six climatic factors using satellite-based products, and investigated the relationships between changes in the lake area and changes in these climatic factors. Some conclusions can be drawn as follows. (1) There was a small numerical difference in the size and change trend of the lake area between multiple scales, and the quasi-monthly lake area had a more subtle, dense and tortuous variation process line. However, whether on the annual scale, seasonal scale or monthly scale, the change of lake area showed a similar significant decline during 1988−2014, and this process can be divided into two sub-stages, namely the first sub-phase with a large area and a slight increase in 1988−1999, and the second sub-phase with a small area and a significant decrease in 2000−2014.
(2) The intra-annual lake area change process showed three types, and it presented the natural bimodal type I in most years in the sub-period of 1988−1999, featuring one low snowmelt runoff recharge peak in March−April and one high rainfall recharge peak in July−September, which is compatible with the natural rhythm of climatic dry−wet conditions; however, as the natural rhythm of the lake changes has been broken by intensive human activities since the late 1990s or early 21st century, the natural bimodal type I has obviously changed into non-natural bimodal type II and unimodal type III, featured by a declining peak in July−September. (3) The climatic dry−wet regime on multi-scales during 1988−2014 in the endorheic Hongjian Lake Basin was generally warming and drying, which is mainly reflected by the increase in temperature (T), arid index (AI) and evaporation (ET 0 , ET a ), and the decrease in the precipitation (Pre) and actual water difference (AWD). There were large differences in the climatic factors at different time scales, especially in the wet and dry seasons. The fluctuation size of climate change at different scales was ranked, in descending order, as the monthly scale, dry season scale, wet season scale, and annual scale. (4) When the lagged effect and cumulative effect were not considered, the correlation coefficient between the lake area and climatic factors was very low, with most values lower than 0.30, and the direction of correlation was not correct; whereas, when gradually considering the lagged effect, the cumulative effect, and the lagged and cumulative combined effect, most correlation coefficients significantly increased to 0.50−0.99, and the direction of the correlation coefficient became coincident with common sense. The correlation analysis identified a lag period of approximately 1−3 years on an annual scale and a lag period of approximately 1−3 months on a monthly scale. This high vulnerability to climate change will mean that lake water management and watershed ecological protection bodies will continue to face high pressures and arduous challenges.