Monitoring of water surface temperature of Eurasian large lakes using MODIS land surface temperature product

In this study, data from MODIS land surface temperature product level 3 (MOD11A2) were used to investigate the spatiotemporal variation of Eurasian lakes water surface temperature (LSWT) from 2001 to 2015, and to examine the most influencing factors of that variation. The temperature of most lakes in the dry climate zone and in the equatorial climatic zone varied from 17 to 31°C and from 23 to 27°C, respectively. LSWTs in the warm temperate and cold climatic zones were in the range of 20 to 27°C and −0.6 and 17°C, respectively. The average day time LSWT in the polar climate zone was −0.71°C in the summer. Lakes in high latitude and in the Tibetan Plateau displayed low LSWT, ranging from −11 to 26°C during the night time. Large spatial variations of diurnal temperature difference (DTD) were observed in lakes across Eurasia. However, variations in DTDs were small in lakes located in high latitude and in tropical rainforest regions. The shallow lakes showed a rapid response of LSWT to solar and atmospheric forcing, while in the large and deep lakes, that response was sluggish. Results of this study demonstrated the applicability of remote sensing and MODIS LST products to capture the spatial–temporal variability of LSWT across continental scales, in particular for the vast wilderness areas and protected environment in high latitude regions of the world. The approach can be used in future studies examining processes and factors controlling large scale variability of LSWT.

LSWTs in the warm temperate and cold climatic zones were in the range of 20 to 27 C and −0.6 and 17 C, respectively. The average day time LSWT in the polar climate zone was −0.71 C in the summer. Lakes in high latitude and in the Tibetan Plateau displayed low LSWT, ranging from −11 to 26 C during the night time. Large spatial variations of diurnal temperature difference (DTD) were observed in lakes across Eurasia. However, variations in DTDs were small in lakes located in high latitude and in tropical rainforest regions. The shallow lakes showed a rapid response of LSWT to solar and atmospheric forcing, while in the large and deep lakes, that response was sluggish. Results of this study demonstrated the applicability of remote sensing and MODIS LST products to capture the spatial-temporal variability of LSWT across continental scales, in particular for the vast wilderness areas and protected environment in high latitude regions of the world. The approach can be used in future studies examining processes and factors controlling large scale variability of LSWT.

K E Y W O R D S
Eurasian lakes, lake surface water temperature, MODIS LST product, spatial-temporal variability 1 | INTRODUCTION Water temperature of inland water bodies is a critical parameter that governs an array of aquatic ecosystem functions that involve physical, chemical, and biological processes (Li, Chen, Zhang, & Recknagel, 2015;Ma et al., 2016). Inland water temperature is regulated by multiple interacting factors, and excessive lake water warming has severe ecological consequences (Kraemer, Mehner, & Adrian, 2017b;Winslow, Read, Hansen, Rose, & Robertson, 2017). Adrian et al. reported that some of the largest lakes across the world have shown notable temperature increases due to climate change (Adrian et al., 2009;Schmidt et al., 2018). For example, the temperature of Lake Baikal (Russia) has increased by approximately 1.2 C since 1946 (Dörnhöfer & Oppelt, 2016). Surface water temperatures in several lakes, including not only large lakes like Superior and Tanganyika but also small lakes like Lower Lake Zurich, are increasing faster than regional air temperatures (Austin & Colman, 2007;O'Reilly et al., 2015;O'Reilly, Alin, Plisnier, Cohen, & McKee, 2003;Schmid & Köster, 2016;Schneider & Hook, 2010). Lake temperature reflects its morphology, watershed conditions and hydrological dynamics, and variably influences the biology of resident aquatic organisms (Edinger, Duttweiler, & Geyer, 1968;Parastatidis, Mitraka, Chrysoulakis, & Abrams, 2017;Song et al., 2016;Wetzel, 2001). The equilibrium temperature and exchange coefficient, which can describe heat exchange processes at the air-water interface, are the key parameters of water temperatures (Edinger et al., 1968;Magee & Wu, 2017;Woolway & Merchant, 2017). Conventional approaches for measuring lake water temperature using in situ sensors and data loggers provide information that is temporally continuous but spatially limited. These limitations therefore hinder widespread application of conventional approaches for capturing temperature variation at large spatial scales (Parastatidis et al., 2017). Yet, synoptic information is needed to capture thermal heterogeneity in large lakes, examine patterns of thermal variation, and explain fundamental biophysical and chemical processes in these water bodies (Hook, Vaughan, Tonooka, & Schladow, 2007;Ke & Song, 2014;. Water temperature of lakes is primarily influenced by their absorption of solar energy in a process affected by an array of physical, chemical and hydro-morphologic properties of lake ecosystems . In this respect, the surface temperature of lake water is highly dynamic, as it changes seasonally and diurnally due to variation in air temperature and the insulation effect of snow/ice Parastatidis et al., 2017;Song et al., 2016;Zhang et al., 2014). The amount of sunlight absorbed by water increases exponentially with the distance it travels through the water column, particularly for radiation wavelengths shorter than 750 nm (Schmidt et al., 2018;Song et al., 2016). The high specific heat of water enables the dissipation of light energy and its accumulation as heat in the water column. However, retention of that energy depends on multiple factors (wind speed, currents and other water movements, watershed geomorphology) influencing its distribution within a lake system, and the change rate between water input and discharge through the tributaries (Gorham, 1964). Thus, the pattern of thermal evolution and stratification influences fundamentally the cycling of physical and chemical components of lakes, which in turn drives primary productivity and decomposition processes (Wetzel, 2001).
Thermal remote sensing methods for monitoring lake surface water temperature (LSWT) circumvent problems of accessibility to lakes in remote areas and so provide a synoptic context for evaluating relationships between landscape features and water thermal characteristics (Kraemer et al., 2017a;Moukomla & Blanken, 2016;Zhong, Notaro, & Vavrus, 2019). Thermal infrared (TIR) remote sensing can reliably map LSWT and its circulation patterns in lakes using various satellite sensors (Ke & Song, 2014). Of those, the Advanced Very High Resolution Radiometer (AVHRR) and Moderate Resolution Imaging Spectroradiometer (MODIS) on board Terra/Aqua have been widely used for monitoring lake surface temperatures (Phillips, Saylor, Kaye, & Gibert, 2016;Zhang et al., 2014). Landsat-Thematic
Although LSWT studies have been conducted at various scales around the world, the distribution pattern of LSWT for Eurasian lakes has yet to be analysed. Investigations about LSWT variation of 246 globally distributed lakes and its climatology have been reported (Layden, Merchant, & MacCallum, 2015). There are ongoing efforts to develop homogenized LSWT from multiple polar orbiting satellites, but these efforts are hampered by the lack of availability for night time (Pareeth, Salmaso, Adrian, & Neteler, 2016). Although there is enough evidence to prove that lakes are warming in global scale, cold and deep lakes often display an amplified response to air temperatures variability (Woolway & Merchant, 2017). Several Investigations indicated LSWT are increasing faster than regional air temperatures by using MODIS LST product (Schneider & Hook, 2010). Torbick, Ziniti, Wu, and Linder (2016) compiled 30 years Landsat satellite data to estimate long-term temperature trend of lakes in the Northeast United States. In the Eurasia, especially Asia, most lakes have no longterm LSWT record and the spatiotemporal patterns are not known.
Two main questions need to be answered: How did LSWT of Eurasian lakes vary, and which factors can contribute to variations in LSWT?
To address these questions, 3,840 composites of MODIS imagery data covering the terrestrial area of Eurasia were processed and analysed. The characteristics of LSWT for 1,098 lakes that are at least 25 km 2 in size were examined to determine spatial associations with climatic, landscape and hydrologic conditions. The objectives of this study were to (a) examine LSWT variation of lakes at the Eurasian

| Criteria of lakes selection
In the pre-satellite remote sensing stage, three types of independent approaches were generally used to obtain a global lake census: (a) a lake-type approach based on the origin of lakes, (b) an extrapolation of known regional censuses, and (c) a climatic approach based on lake distribution in homogeneous temperature and runoff environments (Meybeck, 1995). With respect to the last approach, it links lake distribution with five major climatic biomes, that is, de-glaciated regions (52%), temperate regions (13%), dry and arid regions (25%), deserts (1%), and wet tropics (9%). This approach has clear merits in combining temperature and surface-water runoff (Wetzel, 2001).
With the advent of satellite remote sensing, a global census of lake distribution and thermal characterization has become possible (Downing et al., 2006;Verpoorter, Kutser, Seekell, & Tranvik, 2014).
Using MODIS LST products at a nominal spatial resolution of 1 km 2 , this study extracted lakes and water reservoirs larger than 25 km 2 in Eurasia from the Global Lakes and Reservoirs Database (Lehner & Döll, 2004); this yielded a total of 1,098 water bodies (Figure 1).

Among those, 861 lakes and 237 water reservoirs were in Asia and
Europe, respectively. It should be noted that some lakes or reservoirs greater than 25 km 2 might not have been selected if one of their dimension is smaller than 1-km, the spatial resolution of the MODIS LST product. Such elongated water bodies would likely to be located along shorelines. The selected 1,098 lakes cover 0.63 million km 2 , which is about 44.3% of the world's total lake areas (Lehner & Döll, 2004).

| MODIS LST product pre-processing
Daily MODIS-Terra/Aqua LST (level 2, collection 6) is an option to characterize LSWT. However, in order to track thermal dynamics of water bodies, cloud contamination must be accounted for and it may cause a large proportion of invalid pixels in the MODIS data product.

| LST data analysis
After mosaicking using MRT, the MOD11A2 images were projected onto Albers in GeoTIFF format using nearest neighbour interpolation by using ERDAS Imagine 9.2(ERDAS, 2008). (NIR) and green band, thresholding was used to extract boundaries of water bodies based on the spectral ratio of NIR/Green using Landsat OLI imagery. The shape file for water surface areas greater than 25 km 2 was generated using Landsat-OLI imagery acquired in 2010, by referencing it to MODIS reflectance product (MOD09Q1) in seeking a robust-matching MODIS LST product. To minimize seasonal or inter-annual difference on Eurasian lakes measurement, the algorithm LakeTime, which is suitable for lake mapping at continental and global scales, was used to select Landsat scenes based on seasonally-defined climatic and hydrological variables (Lyons & Sheng, 2018; Sheng F I G U R E 1 Location, across Eurasia, of the 1,098 lakes and reservoirs (blue) selected for the study and the climatic zones where the lakes and reservoirs are situated. Well-known lakes within each climatic zone are shown in red et al., 2016). LakeTime is robust and reliable for delineating lake boundaries using lake-specific NDWI thresholds (Sheng et al., 2016).
LakeTime can not only reduce the impacts of seasonal and interannual variability, but also mitigate cloud contamination effects, and significantly reduce labour costs in the subsequent QA/QC process (Sheng et al., 2016). To avoid possible land-water interface contamination and fluctuation in lake surface areas, a 500-m offshore buffer zone was generated to exclude LST pixels along shoreline zones by using ArcGIS version 10 software (ESRI, 2015;Ke & Song, 2014).
Despite using the 8-day composite product (MOD11A2), invalid pixels could still arise from cloud contamination, so additional procedures were taken to remove questionable pixels in the final MODIS LST product (Deo & Şahin, 2017;Ke & Song, 2014). To do this, unreliable LST pixels over lakes were removed using QC information, ENVI 5.3 and IDL 8.5 (Harris Geospatial Solutions, Inc., 2015), followed by median filtering of time-series LST data stacks (Ke & Song, 2014). In the first step, all pixels having an average LST error of <1 K (i.e., QC = 0, 1, 5, 17, and 21) were kept, while pixels with other QC values were removed through referencing to LST data quality flags stored in the QC file. In the second step, a data cube was generated by stacking 46 extracted LST layers according to the Julian day sequence of the generated product, to which a median filtering algorithm iterated pixel by pixel to eliminate any questionable pixels in the LST images ( Figure 2). In most cases, the 8-day composites consisted of at least 70% valid pixels. Over the 16-year study period, only 97 images (out of 3,840 images or 2.5%) had less than 70% valid pixels. The months with affected images were consistent over the years. Most affected images were from April to June, and from September to mid-December. Furthermore, there were more valid images in the night time MODIS LST product (Ke & Song, 2014).
Finally, to track temporal trends in LSWT for the selected lakes in Eurasia, annual rates of temperature change from MODIS LST data were regressed against the data acquisition year using linear regression (Zhang et al., 2014).
where y represents the MODIS LST, x denotes the time-series of years, e represents the residuals, a is the intercept, and b is the rate of temperature change. The coefficients a and b were determined through least squares fitting. 3 | RESULTS

| Validation of LSWT
One hundred and twenty three LSWT values extracted from MOD11A2 were matched with in situ water temperature measurements at meteorological stations close to Baikal Lake, Vanern Lake, Ladoga Lake and Onega Lake in 2015. The LSWT estimations were spatially averaged over a 3 × 3 1-km pixels window centered on the meteorological stations to achieve spatial representativeness of the measured data. The validation step showed that there is a general underestimation for the four stations, as the bias for all sites is 4.57 K for this site (Figure 3). Figure 3 illustrates the comparison of ground and estimated LSWT, showing that they have a high correlation (R 2 = 0.793). As reported, the overall accuracy of MODIS land surface temperature product was within 1 C on land (Li et al., 2013;Li et al., 2014;Wan, 2014;Wan, Zhang, Zhang, & Li, 2002;Liu et al., 2015). Virdis et al.(2020) reported this level of accuracy can be obtained in an ideal situation, namely when images are acquired with very clear sky conditions (Coll, Wan, & Galve, 2009) and/or when field data are rigorously collected synchronously with satellite passes over permanent monitoring stations specifically built for calibration/ validation activities (Malakar et al., 2018;Sobrino & Skokovi c, 2016). Due to lack of direct in situ LSWT data, our LSWTs compared with MODIS LSWT estimations were land surface temperatures acquired within 200 m close to lake. Although we use 3 × 3 1-km pixels window averaged value to achieve spatial representativeness of the measured data. Variations in the landscape lead to emissivity estimation errors.
It may also be noted that the point scale ground measurement does not exactly correspond to the pixel area retrieved from the satellite.
Considering this limitation, our bias could be acceptable and comparable with the other validation results (Duan et al., 2019;Lu, Zhang, Wang, & Zhou, 2018).

| Patterns of Eurasian lake surface temperatures
As shown in Figure   Note: Area, lake area at the surface (km 2 ); DT, day time temperature ( C); DTD, diurnal temperature difference ( C), Latitude (degree), Water storage (km 3 ); NT, night time temperature ( C); Water depth, maximum water depth (m).
F I G U R E 5 Inter-annual temperature variations of 10 selected well-known large lakes from Asia. Please note that LSWT is plotted on different scales for each lake area may have more influence on the DTR. This issue is worth further study.

| Inter-annual variation in lake surface temperatures
Following the IPCC AR.5 report (IPCC, 2014) and the recommendation "to track global warming hotspots with reference to global climatic zones and corresponding lake distributions," a group of 10 wellknown lakes (largest lakes within a climatic zone) from Asia and 8 well-known lakes from Europe ( Figure 1) were selected, and analysed with respect to inter-annual variation in LSWT using MODIS LST measurements from 2001 to 2015.

| Well-known lakes from Asia
The inter-annual variations in LSWT for the 10 well-known lakes from  (Figure 5a). On the Mongolia Plateau, a slight but insignificant rise in LSWT was apparent for Lake Khar Us at night time, whereas Lake Uvs showed a slight temperature decrease, albeit not significant, over the study period (Figure 5e and Figure 5b). Likewise, for Lake Balkhash and Khanka, a slight decrease at day time and a slight increase at night time were observed, but these trends were not significant (Figure 5c and Figure 5f). Situated in the Pamir, the alpine Lake Issyk underwent a rather weak temperature decline over time, but this was not a significant trend (Figure 5d). It is worth noting that Lake Taymyr (Figure 5h), did show a strong and obvious increase in temperature during the study period (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) (Figure 5j).

| Well-known lakes from Europe
The inter-annual variations of eight selected well-known lakes across F I G U R E 6 Inter-annual temperature variations of 8 selected well-known large lakes from Europe. Please note that LSWT is plotted on different scales for each lake 4 | DISCUSSION

| Ice-free line for Eurasian lakes
Significant variability (S.D. = 187.9 days) exist among Eurasian lakes with respect to the duration of the ice-free period, which can last between 24 and 365 days (Figure 7). At the continental scale, ice-free duration had a clear latitudinal pattern that matched the air temperature distribution. Lakes located in the boreal and tundra biomes, along with some lakes located at high elevation, exhibited short durations of ice-free period (24 to 168 days). Conversely, with the exception of those located at high elevation, lakes in the tropics were generally icefree during the whole year. In general, lakes in the Tibetan Plateau had shorter ice-free durations compared to lakes located at lower elevation but at the same latitude. Everything else being the same, lakes with a large water volume generally had a longer ice-free duration than smaller lakes.

|
Regulating factors of lake surface temperature 4.2.1 | Climatic zone and air temperature LSWT varies according to climate zones (Figure 8). The temperature of most lakes in the dry climate zone varied between 17 and 31 C, while for most lakes in the equatorial climatic zone LSWT was between 23 and 27 C. In the warm temperate climatic zone, LSWT was in the range of 20 to 27 C. In the cold temperate zone, for most lakes (620 lakes or 56.5% of Eurasian lakes) LSWT was between −0.6 and 17 C. The average summer day time LSWT was −0.71 C in the polar climate zone.
The water temperature of lakes on the Tibetan Plateau is generally lower than in lakes at the same latitude but in other climatic regions. As shown in Table 1, the elevation of the three Tibetan lakes selected was >4,000 m, while the elevation of the three selected lakes in the Indo-Pakistani region was below 600 m (although at the same latitude). So, lakes elevation is a critical factor explaining difference in LSWT. Based on our results, for every 100 m increase in elevation, the LSWT drops by about 0.6 C. The higher the altitude, the lower the lake water surface temperature.
All of the large water bodies in the South Asian subcontinent are man-made water reservoirs (Figure 8). In the Deccan Plateau (India), the altitude of the three selected reservoirs exhibiting highest and lowest LSWT was between 425 and 519 m. Land use type in the watersheds surrounding these reservoirs was generally similar (paddy field). Water depth was also similar (Table 1). So, the key difference among these three reservoirs was their climatic zones location. Therefore, the main driving factor of LSWT in the South Asian subcontinent is climate zone, rather than lake morphometric characteristics.
Regression models were generated between in-situ air temperature and lakes LSWT derived from MODIS to examine relationship between these variables (Figure 9). Positive relationships were found

| Lake area versus lake depth
Area, depth and volume of lakes have strong impact on LSWT (Figure 3, Table 1). Of the lakes investigated, Lake Tonlesap (Cambodia) exhibited the highest temperature (max: 26.87 C); it is the shallowest of the lakes (around 12 m deep) and is located at the lowest elevation (4 m) ( Table 2).
A previous study has suggested that the annual heat budget of lakes is strongly dependent on mean depth because a deeper water column corresponds to a larger heat storage capacity (Gorham, 1964). That study further indicated a strong relationship between heat budget and lake water volume (semi-logarithmic scale) (Gorham, 1964). For deep and large lakes, more energy is needed to raise lake water temperature by one unit. Thus, less temporal variation in LSWT was observed with the very large lakes, for example, Baikal, Onega, and Issyk (Figure 5 and

| Ice cover trend and salinity
On a continental scale, lakes at high latitudes and high altitudes tend to have low LSWT and short ice-free duration (Figure 7 and Table 3).
Lakes Niagame and Kungasalakh, located near the Arctic Circle, exhibited much shorter ice-free time (average: 56 and 64 days, respectively) compared to Lake Eloygytgyn and Yambuto in eastern Siberia (128 and 136 days, respectively). Since these lakes are found at similar elevations and are located in the same climatic zone, water chemistry could be a reason for the difference in ice-free duration. It is well known that increased salinity of lake waters can increase ice formation temperature, thus an extension in ice-free duration should be expected for lakes with high salinity (Wang & Dou, 1998). The merit of that contention cannot be determined by the data presented in this study, but it does underscore the need to further investigate how salinity could affect the dynamics of ice formation, the duration of the ice-free period in lakes, and ultimately the response of water temperature in these water bodies to a warming global climate.

| CONCLUSIONS
The study demonstrated the utility of MODIS LST products to monitor the spatial, diurnal, seasonal, and annual variations of lakes water F I G U R E 9 Relationships between lake water surface temperature (LSWT) and air temperature during the day time (a) and night time (b)

T A B L E 3
The top ten lakes (from a set of 1,098 lakes) with the lowest water surface temperature across Eurasia Note: ASL, elevation above sea level (m); DT, day time temperature ( C); DTD, diurnal temperature difference ( C); NT, night time temperature ( C).
surface temperature (LSWT) across Eurasia where long-term in-situ measurements are very limited. In particular, this is an effective approach for monitoring LSWT in some of the wild, inaccessible and