**3. Methods**

### *3.1. Landsat Data and Lake Surface Area Extraction*

Long-term changes in the surface area of Tuosu, Keluke, and Gahai Lakes were extracted from remote sensing images. Since 1984, Landsat satellites have acquired highresolution Earth observation images, which are widely used for feature identification. In this study, the Landsat 5 Thematic Mapper (TM) and Landsat 8 Operational Land Imager (OLI), which are the longest time-series data currently available, provided observation data for different periods. The Universal Transverse Mercator (UTM) and World Geodetic System 1984 were used as geocoordinate references to construct the temporal and spatial sequence of lake changes. The spatial resolution of the extracted lake surface area was 30 m, which is sufficient for studying lakes measuring several tens to hundreds of square kilometers. Although the spatial resolution is moderate, its effect on the results is very limited [22–25]. All Landsat data were downloaded from the United States Geological Survey (http://glovis.usgs.gov/, accessed on 10 January 2021), Geospatial Data Cloud (https://geocloud.cgs.gov.cn/, accessed on 10 January 2021), and Chinese Academy of Sciences (http://www.gscloud.cn/, accessed on 10 January 2021). The necessary image preprocessing steps, such as radiation calibration and atmospheric correction, were performed using ENVI 5.3 software.

The repeat coverage of Landsat 5 and Landsat 8 is 18 and 16 days, respectively, providing 1–2 images per month. Considering the influence of seasonal variations, monthly variations of the Tuosu Lake surface area in 2003 and 2020 were also extracted (Figure 2), and it showed that the seasonal variation had only a slight disturbance to the annual trend of lake surface area (see Section 4.1 for details). To observe the lake surface area in more detail, this study selected 20 remote sensing images of the Delingha region from 2000 to 2020 and extracted changes in the surface areas in Tuosu, Keluke, and Gahai Lakes. As optical images are significantly affected by weather conditions, particularly clouds, satellite images were typically only selected from July to August. This is because the summer in Qaidam Basin is dry and sunny with minimal clouds. If no images were available, we substituted images from adjacent months. The data sources used in the study are listed in Table 1. It should be noted that there was a lack of images available in 2012.

**Figure 2.** Changes in the surface area of Tuosu Lake in 2003 and 2020. Black error bar represents the uncertainty of the extracted lake surface area.

**Table 1.** Sources of remote sensing data used in this study.


Uncertainty in the extracted lake surface area mainly originates from the positioning accuracy and indistinguishable mixed pixels in the image [26]. According to previous studies, the registration error is 6 m for TM images [27] and 5 m for OLI images [28]. Considering the registration error, the positioning accuracy in the image, and the clarity of the lake boundary, the uncertainty in the lake surface area was estimated by Equation (1) [29] and is shown as error bars in Figure 3a. Here, *EA* is the uncertainty in the extracted lake surface area; *l* is the length of the lake boundary; *LREyear* is the resolution error of Landsat images in different years, which should be half the resolution of the image pixel; and *Eco* is the registration error in the image:

$$E\_A = l \times \sqrt{LRE\_{year}^2 + E\_{co}^{-2}} \tag{1}$$

Because of a lack of detailed bathymetric maps of the lakes, the change in lake storage was estimated based on the lake surface area and the slope of the lakeshore zone. The change in water storage between the two stages can be approximated as a frustum, and its volume can be estimated by the following equation [30]:

$$
\Delta V = \frac{1}{3} \times \left( S\_1 + S\_2 + \sqrt{S\_1 \times S\_2} \right) \times H \tag{2}
$$

where *S*1 and *S*2 are the lake surface area in the two stages; Δ*V* is the change in lake storage; and *H* is the lake water level interval. The lakes in the study area are located at the end of the alluvial lacustrine plain, with gentle topography and a lakeshore slope of less than 5%, where the slope *I = H/L*, and *L* is the horizontal distance. For *H*, the shapes of Tuosu, Keluk, and Gahai Lakes are approximately equilateral triangles, whose side lengths can be calculated according to the lake surface area in different periods:

$$H = \frac{2}{3} \times I \times \left(\sqrt{S\_2} - \sqrt{S\_1}\right) \tag{3}$$

here, slope values of 4, 3, and 2 were used to calculate the change in lake water storage, and the calculation results are shown in Figure 3b.

**Figure 3.** (**a**) Plot depicting the change of lake surface area in the study region extracted from Landsat satellite images. Triangles (before 2000) represent data in the literature [31,32]. y, y1, and y2 represent the least square fitting line of the surface area in Gahai Lake after 2000, Tuosu Lake from 2003 to 2011, and after 2011, respectively. Black error bar represents the uncertainty of the extracted lake surface area. (**b**) Calculated change in lake storage (Δ*V*). The solid line, dashed lines in dark gray and light gray are the calculation results when the slope (*I*) is 3, 2, and 4, respectively.

### *3.2. Hydrological and Meteorological Monitoring Data*

The meteorological monitoring network in the alpine region is sparse; however, the single weather station is generally considered to represent the climatic conditions of the basin. The distance between Delingha station (97◦22E, 37◦22N) and Tuosu, Keluke, and Gahai Lakes is only 42, 41, and 27 km, respectively. The terrain between the weather station and the lakes is flat, with no obstruction from mountains, and the climate conditions exhibit minimal spatial variability. Therefore, data from Delingha station were used to represent the climatic conditions near the lakes. Meteorological data at Delingha station were downloaded from the National Meteorological Data Center (http://data.cma.cn/, accessed on 9 August 2021); this included daily data of relative humidity, minimum relative humidity, mean temperature, maximum temperature, minimum temperature, mean wind speed, minimum wind speed, and annual precipitation. Delingha station lacks solar radiation data; therefore, these data were replaced by monitoring data from the nearest station in Golmud (94◦54E, 36◦25N).

Runoff data from Delingha station (1957–2018) and Zelinggou station (1957–1983) were compiled from the *Hydrological Yearbook*. Lake water temperature data for Tuosu, Keluke, and Gahai Lakes were obtained from a dataset of daily lake surface temperature over the Tibetan Plateau (1978–2017) compiled by the National Tibetan Plateau Data Center (http://data.tpdc.ac.cn, accessed on 11 August 2021) [33], which uses the improved lake water temperature model (air2water) to simulate the annual surface temperature every day. Considering the different available periods of sequences, the calculation interval for analyzing lake evaporation was 1984–2016 in this study.

### *3.3. Penman–Monteith Model*

Accurate quantification of lake evaporation is essential for determining the water budget of a lake, especially for endorheic lakes, where evaporation is the most important

output term of the lake water budget [34–36]. Based on the energy conservation formula, Penman [37] first proposed the Penman formula for calculating evapotranspiration, from which many water surface evaporation models have been developed [38,39]. After comparing 14 evaporation models with baseline Bowen ratio energy budget measurements, Rosenberry [40] identified that the De BruinKeijman, Priestley Taylor, and Penman models provided the best estimates of water surface evaporation. McJannet then compared eddy covariance measurements with the De BruinKeijman, Priestley Taylor, and Penman– Monteith models, [15,16] and found that the Penman–Monteith model was most suitable for water surface evaporation because it considers both the vapor pressure gradient and wind speed, and produced estimates of total evaporation that varied from the actual measurements by less than 1%. In addition, the Penman–Monteith model has been improved by proposing the general application of wind functions, making it applicable for calculating evaporation for water bodies ranging from tens to hundreds of kilometers. This improved Penman–Monteith model has been applied in several subsequent studies [17–19].

In this study, we used the improved Penman–Monteith model to calculate lake surface evaporation (see McJannet [15,16] for details). Uncertainty in the evaporation value mainly derived from data measurement and parameter calculation. The measurement error comes from the measurement of wind speed, temperature, humidity, and accumulated solar radiation. These errors are inevitable but have little impact on the final calculation results.

### *3.4. Sampling and Isotope Measurements*

Hydrogen (2H) and oxygen (18O) isotopes are widely used to study hydrological cycles [41,42] and qualitatively identify water sources and trace groundwater runoff processes [11,13]. In this study, 47 water samples were collected from rivers, phreatic groundwater, confined groundwater, springs, and lakes (Figure 1). River water samples (R07, R09, and R12) were collected from the upper reaches of Bayin River, close to where the river exits the mountains, mainly from glacial meltwater; thus, they are minimally influenced by precipitation and groundwater. These three samples represent the hydrogen (2H) and oxygen (18O) isotope characteristics of glacial meltwater. In this study, only one lake water sample was collected from Keluke Lake, which was analyzed together with 12 lake water samples collected from Yang Lake [11], which are discussed in Section 4.4.

Stable isotopes of hydrogen and oxygen (2H, 3H, 18O) were measured at the State Key Laboratory of Hydrology, Water Resources, and Hydraulic Engineering, Hohai University. 18O/16O and 2H/1H ratios were measured using a MAT253 mass spectrometer. The isotope ratio '*δ*' was expressed as follows:

$$
\delta\_{sample}(\%) = \left( R\_{sample} / R\_{standard} - 1 \right) \times 1000 \tag{4}
$$

where *Rsample* and *Rstandard* are the isotope ratios (18O/16O, 2H/1H) of the sample and standard, respectively, and the international standard is the δ2H and δ18O of Vienna mean seawater. The measurement errors of δ18O and δ2H were ±0.1 ‰ and ±1‰, respectively. The tritium content of samples was measured by a liquid scintillation meter (TRI-CARB 3170 TR/SL) with a detection limit of 0.2 TU and precision of >0.8 TU. The measurement results for hydrogen (2H, 3H) and oxygen (18O) isotopes are shown in Table 2.



Note: PGW, CGW, RW, SW, and LW represent phreatic water, confined water, river water, spring water, and lake water, respectively.
