Lake Surface Water Temperature Derived from 35 Years of AVHRR Sensor Data for European Lakes

Abstract

Lake surface water temperature (LSWT) is an important parameter with which to assess aquatic ecosystems and to study the lake’s response to climate change. The AVHRR archive of the University of Bern offers great potential to derive consistent LSWT data suited for the study of climate change and lake dynamics. To derive such a dataset, challenges such as orbit drift correction, non-water pixel detection, and homogenization had to be solved. The result is a dataset covering over 3.5 decades of spatial LSWT data for 26 European lakes. The validation against in-situ temperature data at 19 locations showed an uncertainty between ±0.8 K and ±2.0 K (standard deviation), depending on locations of the lakes. The long-term robustness of the dataset was confirmed by comparing in-situ and satellite derived temperature trends, which showed no significant difference. The final trend analysis showed significant LSWT warming trends at all locations (0.2 K/decade to 0.8 K/decade). A gradient of increasing trends from south-west to north-east of Europe was revealed. The strong intra-annual variability of trends indicates that single seasonal trends do not well represent the response of a lake to climate change, e.g., autumn trends are dominant in the north of Europe, whereas winter trends are dominant in the south. Intra-lake variability of trends indicates that trends at single in-situ stations do not necessarily represent the lake’s response. The LSWT dataset generated for this study gives some new and interesting insights into the response of European lakes to climate change during the last 36 years (1981–2016).

1. Introduction

Water temperature is an important parameter within aquatic ecosystems. The metabolic rate of organisms and the activity of bacteria and algae is governed by water temperature [1,2,3,4]. Changing water temperature can alter the balance and the species composition of an aquatic ecosystem and even alter its trophic state [5,6,7]. Lake dynamic processes such as thermal stratification and mixing dynamics are firmly related to water temperatures [8,9]. The energy balance of large lakes is dominated by heat exchange processes at the surface [8]. Lake surface water temperature (LSWT) governs energy and mass exchange between the water body and the troposphere. It is also an important parameter in mesoscale meteorology and local climate models [10,11,12]. LSWT can be seen as an indicator of regional climate change [13], and several recent publications analyzed climate signals within lakes by investigating LSWT trends [14,15,16,17]. Recognizing its importance, LSWT was also added to the list of essential climate variables (ECV) by the Global Climate Observing System (GCOS), which needs to be monitored [18].

Traditionally, measurement of water temperature is done with in-situ thermometers, either manually or with automatic recording stations. Depending on the objectives of organizations or institutes, measurements are taken at various point-locations, at various intervals, and at various depths. This leads to a high heterogeneity amongst available lake water temperatures data, preventing the composition of consistent time series to be analyzed in the frame of climate change studies. Also, for lake dynamic studies, the 1D point measurements from automatic in-situ stations are not well suited; data with spatial coverage would be more valuable.

Remotely sensed earth observation data from satellites has the potential to overcome these limitations. For the study of lake dynamic processes and climate change, the spatial and temporal homogeneity of a dataset is of crucial importance, whereas the accuracy of an individual point-measurement is secondary. Satellite observation does not provide highly accurate temperature measurement for a single point in space. It provides a grid of temperature measurements, each averaged over a certain portion (depending on the resolution) of the lake’s surface. It is also supposed to have higher uncertainties due to atmospheric disturbance. However, the potential of temporal and spatial coverage is unique, and therefore remotely sensed surface water temperature from satellites is a highly relevant topic in the scientific community.

Having the WMO definition of a climate period in mind, homogeneous time series covering more than 30 years would be necessary to be relevant for the study of climate change. The Advanced Very High Resolution Radiometer (AVHRR), mounted on various satellite platforms since 1981, offers the possibility to generate such a data set. In 1981, the first AVHRR sensor recording in two TIR bands was launched with the NOAA-7 platform. Since then, the AVHRR sensor was mounted on the platforms of the NOAA series (until NOAA-19), and after on the last MetOp series, operated by EUMETSAT. This leads to an archive of AVHRR data starting in June 1981 until today. This time series, which covers over 35 years, is unique in satellite-based remote sensing, and therefore, despite the moderate resolution and the noise rate, the almost antique AVHRR sensor is still of interest today.

For this study, we generated a LSWT dataset of European lakes covering 36 years from 1981 to 2016 considering AVHRR local area coverage (LAC) data, archived at University of Bern, with a spatial resolution of 1.1 km in nadir. This data allows one to gain a unique insight into the response of European lakes to climate change in the last 3.5 decades.

The next chapter gives an overview of data and study sites considering 26 lakes distributed in Europe; in Section 3, the processing method is described in detail, which includes also post-processing based on different quality indicators and a correction of the diurnal temperature cycle to minimize the influence of orbit drift. Section 4 describes the results of this study, including validation with in-situ measurements and, finally, shows the trends for European lakes. The discussion in Section 5 compares the findings of LSWT trends with air temperatures of different periods and shows the influence of the measurement site on the trend. Section 6 concludes the paper.

2. Data and Study Sites

2.1. The Study Site

Lakes in and around Europe were selected to have a good distribution from east to west and from north to south (Figure 1). Only lakes exceeding a certain size were retained, because multiple adjacent water-only pixels must be available to perform accurate LSWT retrieval. In the case of the AVHRR sensor with a pixel size of roughly 1 km, this means that a lake should contain patches with open waters exceeding 3 × 3 km. For validation purposes, lakes with suitable in-situ temperature measurements also had to be among the selected lakes, in which ‘suitable’ means that the measurement should be representative of the open surface waters of the lake in terms of location and measurement depth.

With these criterions, a preselection adapted from Politi et al. [19] was used first, followed by a final manual selection. This led to a selection of 26 lakes (

Table 1. List of Lakes where LSWT was retrieved and trend analysis performed. The acronyms are the same as those used in Figure 1.

Lake Name	Country	Acr.	Location		Alt. [m]	Area [km2]	Depth Ø [m]	In-Situ Data
Lake Balaton	HU	Bala	46.75°N	17.43°E	126	578	3.3	Yes
Lake Bolsena	IT	Bols	42.61°N	11.93°E	368	114	81.0	No
Lake Constance	DE,CH,AUT	Cons	47.60°N	9.45°E	431	472	101.0	Yes
Djorf Torba	DZA	Djor	31.53°N	2.75°W	716	77	6.1	No
Lake Garda	IT	Gard	45.52°N	10.67°E	85	368	133.0	No
Lake Geneva	FR, CH	Gene	46.41°N	6.41°E	466	580	152.7	Yes
Inarijärvi	FIN	Inar	69.04°N	27.76°E	126	1′110	14.0	Yes
Kiev Reservoir	UA	Kiye	50.80°N	30.47°E	93	420	4.0	No
Lappajärvi	FIN	Lapp	63.13°N	23.71°E	85	141	7.0	Yes
Lokka Reservoir	FIN	Loka	67.94°N	27.56°E	259	488	5.0	No
Müritz	DE	Muri	53.43°N	12.67°E	63	110	6.5	Yes
Lough Neagh	NIR	Neag	54.62°N	6.41°W	58	380	9.0	No
Lake Ohrid	MK,ALB	Ohri	41.02°N	20.72°E	700	347	155.0	No
Lake Onega	RU	Oneg	61.63°N	35.53°E	56	9′608	30.0	No
Oulujärvi	FIN	Oulu	64.29°N	27.22°E	128	890	7.6	Yes
Lake Päijiänne	FIN	Paij	61.52°N	25.42°E	95	1′080	16.0	Yes
Lake Peipus	EST,RU	Peip	58.53°N	27.50°E	29	3′553	7.1	Yes
Lake Pielinen	FIN	Piel	63.24°N	29.64°E	113	870	9.0	Yes
Säkylän Pyhäjärvi	FIN	Pyhj	61.00°N	22.28°E	52	150	5.5	Yes
Pyhäselkä	FIN	Oriv	62.60°N	29.77°E	76	361	8.7	Yes
Sanguinet	FR	Sang	44.49°N	1.15°W	24	57	5.0	No
Thingvallavatn	IS	Thin	64.16°N	21.13°W	105	83	34.0	No
Lake Trichonida	GR	Tric	38.57°N	21.57°E	40	95	30.0	Yes
Lake Vänern	SWE	Vane	59.07°N	13.60°E	45	5′650	27.0	No
Lake Vättern	SWE	Vatt	58.33°N	14.50°E	91	1′890	40.0	No
Lake Yalpuh	UA	Yalp	45.42°N	28.63°E	10	225	2.0	No

), with latitudes reaching from 31.53°N (Djorf Torba reservoir in Algeria) to 69.04°N (Inarijärvi in Finnland), and longitudes from 21.13°W (Thingvallavatn in Iceland) to 61.90°E (Lake Onega in the Russian Federation). The exact locations, where the LSWT data was resampled for the final time series, are shown in Table The sampling points were preferentially placed in the middle of large undisturbed open water surfaces, as the priority was to get the best possible conditions for an accurate measurement. Whereas, for the validation with in-situ data, sometimes different locations within the lakes were chosen (

Table 2. List of in-situ measuring locations from which data was available for this study. Most of the measurements are either hourly or daily averaged temperatures from automatic stations. In the case of Müritz, Lake Trichonida, and the older measurements from Lake Balaton, the measurements are sporadically and manually taken at specific locations on the lake. The average measurement density varies from 12 measurements per year at Lake Trichonida (monthly measurements) to 8760 measurements per year at Bregenz, Lake Constance (hourly measurements).

Lake	Country	Sta. Name	Type/Freq.	Coverage	N/Yr	Location
Lake Balaton	HU	Keszthely old	instant	1981–2003	21	46.74°N 17.29°E
		Keszthely	hourly	2003–2013	4′043	46.76°N 17.27°E
		Siofok	daily	1992–2012	365	46.90°N 18.05°E
Lake Constance	DE, CH AUT	Bregenz old	daily	1989–1996	365	47.51°N 9.74°E
		Bregenz new	hourly	1997–2009	8′760	47.51°N 9.74°E
		Lindau	hourly	2005–2014	8′261	47.54°N 9.68°E
Lake Geneva	FR, CH	Buchillon	hourly	2000–2011	8′147	46.46°N 6.40°E
		INRA	daily	1991–2011	360	46.38°N 6.47°E
Lake Inari	FIN	Inari	daily	1980–2014	153	68.90°N 28.10°E
Lappajärvi	FIN	Halkosaari	daily	1980–2014	167	63.26°N 23.63°E
Müritz	DE	Müritz	instant	1980–2013	18	53.42°N 12.70°E
Pyhäjärvi	FIN	Säkylän	daily	1996–2014	205	61.11°N 22.17°E
Oulujärvi	FIN	Oulujärvi	daily	1980–2014	156	64.40°N 27.17°E
Lake Paijiänne	FIN	Päijätsalo	daily	2003–2014	200	61.47°N 25.59°E
		Kalkkinen	daily	1980–1999	214	61.29°N 25.57°E
Lake Peipsi	EST, RU	Mehikoorma	daily	1985–2011	267	58.24°N
Lake Pielinen	FIN	Nurmes	daily	1980–2014	163	63.54°N
Pyhäselkä	FIN	Joensuu	daily	1980–2014	194	62.60°N
Lake Trichonida	GRE	Trichonida	instant	2004–2006	12	38.57°N

) to be as close as possible to the in-situ measurement location. These validation points were therefore often close to the shores, which had the effect of decreasing the availability and reliability of the satellite data while increasing the comparability with the in-situ measurements.

2.2. The Data

The AVHRR data archive from the University of Bern [20] was used to retrieve LSWT time series. The archive contains around 130,000 scenes from the European subset (Figure 1) and includes data from all NOAA and MetOp satellites with an AVHRR-2 or AVHRR-3 sensor on board (Figure 2). The data in the archive was preprocessed (see [20] for details) and included a probability cloud mask after Musial et al. [21] for each scene. The data of NOAA-15 was omitted, because the data from this platform showed significantly increased amounts of erroneous pixels and scanlines. The data from MetOp-1 was also omitted due to inconsistencies during the pre-processing.

The in-situ data available for this study is shown in

Table 3. The six quality levels indicate which additional tests are passed. The levels are cumulative (e.g., Q4 means Q3 AND VZA < 55). The average data availability (compared to Q1) and the average standard deviation during the validation are shown for the different quality levels.

Quality Level		Perc. N (%)	Stddev
Q1	valid pixels *	100.0%	1.399
Q2	SSDEV < 3	99.9%	1.395
Q3	no Sun glint	95.9%	1.390
Q4	VZA < 55	89.7%	1.362
Q5	VZA < 45	83.7%	1.344
Q6	feature detection passed **	66.4%	1.333

* Cloud and land mask applied, and physical inconsistency removed; ** passed if the detection of the feature is confirmed (additional georeferentiation test).

are, mainly, two types of measurements: homogeneous averaged measurements on a daily or hourly basis, and sporadic instant measurements at a defined point in time. For the measurements at lake Trichonida and Lake Balaton (Keszthely old), only the date (day) without exact timestamp (hh:mm) was available; therefore 12:00 local time was assumed. All the measurements are representing bulk temperatures close to the surface (taken at 0.1–1.0 m). Further, most of the northern lakes are ice covered during a certain amount of time per year. For these lakes, the amount of data per year (N/Yr. in Table 2) is considerably below of the expected 365 measurements (daily averaged data) per year.

The in-situ data were filtered for outliers and physically inconsistent values, as well as checked for irregularities such as systematic temperature changes due to modification of the equipment or the location. In the case of the INRA (Institute Nationale de la Recherche Agronomique) measurement at Lake Geneva, a systematic temperature decrease of roughly 0.4 °C was detected after autumn 2007, probably due to a dislocation of the station. Therefore, only the data from 1991 to 2007 from the INRA station was used.

It is assumed that the best validation performance will result from comparison with the hourly averaged and instant measurements, whereas for the daily measurements, the error between an individual satellite measurement and the daily average temperature can theoretically be as high as the amplitude of the diurnal temperature cycle.

3. Method

3.1. LSWT Processing Chain

Figure 3 shows a scheme of the main LSWT retrieval and processing steps performed for this study. For accurate LSWT retrievals, also the pre- and post-processing of the satellite data is very important, e.g., only reliable clear sky observations with accurate geocoding should be retained to ensure that the resulting temperatures correspond to the features we actually investigate. Conservative cloud and land masking are responsible for the elimination of non-water pixels. Parameters like view angles and spatial consistency of LSWT are additional indicators for the quality and reliability of the measurements. All this information is retained in a quality mask associated with the final LSWT product.

The pre-processing of the AVHRR data used for this study is described in Hüsler et al. [20] and provides georeferenced AVHRR scenes, as well as the associated cloud mask after Musial et al. [21]. The coefficients used for this study are time-independent and lake-specific coefficients, empirically matched to radiative transfer model results as described in Lieberherr et al. [22]. The details of the final LSWT retrieval is described within the following subchapters.

3.1.1. Georeferencing Test

Even on well geocoded AVHRR scenes, the accuracy of the geolocation can vary within the same scene. Thus, an additional georeferencing test was added to the pre-processing described in Hüsler et al. [20]. During the extractions of the sub-sets from the scene, an algorithm verifies and adjusts the local geocoding for each target location. The test is basically an edge detection and matching algorithm. First, an edge detection with a Laplace filter is applied to the subset of the scene, and to the subset of a reference image. Then, the correlation for each possible positioning of the AVHRR subset within the reference subset is computed, resulting in a correlation map. In the case of a correct georeferenced image, a correlation peak is found in the center of the correlation map. Otherwise, depending on the position of the correlation peak, a latitude and/or longitude shift is performed before the subset of the lake is extracted. To determine if the correlation maximum corresponds to the correct location, it is tested if the correlation peak is unique and of significant amplitude. The geolocation check is performed with AVHRR Band 2 data for daytime scenes, and with AVHRR Band 4 data for night-time scenes.

3.1.2. The LSWT Retrieval

A state-of-the-art split-window method was used to retrieve surface water temperature for this study. The idea is to make use of the different absorption characteristics of water vapor within differential spectral windows of the thermal infrared [23]. Two algorithms are currently used in various applications: the Multichannel Sea Surface Temperature (MCSST) algorithm [24] assumes a linear relationship between the atmospheric transmittance of the signal and the difference in absorption of the two channels. In contrast, the Non-Linear Sea Surface Temperature (NLSST) algorithm introduced in 1990 for NOAA SST products takes the nonlinear nature of atmospheric absorption into account, which increases the accuracy especially for daytime measurements [25,26]. Recent publications show no significant improvement when using the NLSST instead of the MCSST for LSWT retrieval [27,28,29]. Therefore, here the MCSST approach described by Equation (1) with empirically determined split window coefficients is used to derive LSWT from AVHRR data.

$$\mathrm{LSWT}=a_0+a_{1}B{T}_4+a_2\left(B{T}_4-B{T}_5\right)+a_3(1-\sec \theta )\left(B{T}_4-B{T}_5\right)$$

(1)

• $a_{0,1,2,3}$ Split-window coefficients
• $B{T}_{4,5}$ Brightness temperature from AVHRR band 4 and 5
• $\theta$ Satellite view zenith angle

The coefficients where derived for each location (lake specific) and satellite, based on the findings of Lieberherr et al. [22]. It was assumed that the use of matchup data from a 2-year period is sufficient to obtain the infinite temporal validity of the coefficients. The used matchup data was generated with a radiative transfer model, RTTOV-11 [30], using atmospheric profiles of water vapor and temperature from the ERA-INTERIM [31] dataset.

3.1.3. Non-Water Pixels Elimination

Cloud and land masking is crucial for satellite-based temperature retrievals of surface lake water. Clouds are opaque to thermal infrared radiation, and even a slightly cloud-contaminated signal (sub-pixel sized clouds) can result in a considerably colder measurement of the temperature. Also, for land contaminated pixels, the signal can be biased, as the land surface temperature sometimes differs by several degrees from the water temperature. Therefore, strict cloud and land masking are necessary to eliminate affected pixels. On the other hand, this reduces the amount of available data to be used for further processing. An ideal cloud and land masking would eliminate all affected pixels without any false positives. As, of course, this is not feasible, a good trade-off tailored to the application must be found. For the cloud detection within this study, the probability cloud mask developed by [21] is used. Additionally, in order to improve the performance over inland-waterbodies and during night-time, some additional tests are performed. These additional tests are described in [32]. The masking of land pixels is done by applying a land-water mask (LWM).

The LWM is derived from freely available Open Street Map (OSM) data, which was found to be more accurate in geolocation than several other land classification datasets. The accuracy of OSM data was the subject of various studies, e.g., [33] reported an accuracy of ±3 m (SD) in an urban area in Germany, and [34] found an RMSE of 15 m in a remote rural area in Portugal. Further, the reliability of OSM data is considered good, since it is based on in-situ GPS measurements and is verified and reviewed by the OSM community. For easily accessible and frequently visited objects such as Lake Constance, the community is very active, and the data is highly reliable. When rasterizing the OSM-shapefile to the AVHRR pixel grid, the LWM is extended by 500 m (~1/2 pixel) into the lake to exclude lake pixels where the signal is potentially land contaminated. The low number of lakes (26) allowed us to check the thematic and geographic accuracy manually during the generation of the LWM raster files from the OSM water shapes.

3.1.4. Quality Levels

The reliability of each measurement (pixel) is expressed by a quality score reaching from Q0 to Q6; Table 3 recapitulates the six quality levels. The tests for the quality levels are cumulative, meaning that for each level, the tests of the previous quality levels have been passed too. Pixels with physically inconsistent values (LSWT > 35 °C or LSWT < −5 °C), as well as isolated pixels, are assigned to Q0; all other clear-sky pixels and water pixels are classified at least QPixels with a low spatial standard deviation (similarity to surrounding pixels) get QPixels outside the potential sun glint zone are assigned QQ4 and Q5 are achieved if the view zenith angle (VZA) is lower than 55° or 45°, respectively. Q6 is assigned if the georeferentiation is considered highly reliable. This means that during the georeferencing test, the feature (lake) must be detected with certainty. The results showed that at Q5 and Q6 (VZA < 45 and the enhanced georeferentiation test), the available pixels decrease strongly without considerably reducing the uncertainty. Finally, Q4 was used for the trend analysis during this study, as it was considered to be the best tradeoff between accuracy and data availability.

3.2. Time Series Generation

Once all the AVHRR scenes are processed, the 1D temperature timeseries for each specific target point are extracted from the 2D LSWT product. For each target point, pixels within a radius of 3 km are extracted. All invalid pixels and all pixels with a quality below the aimed quality level are removed. If at least 3 pixels remain within the 3 km radius around the target, the average of the valid pixels is retained. The thus-sampled temperature value is then integrated into the raw time series with the timestamp of the satellite scene accorded.

Before using this raw time series for validation and for trend analysis, an outlier detection adapted for such time-series is applied. The main problem with basic outlier detection on temperature time series is that the yearly and the diurnal temperature cycle already induces a high amount of variance in the data. For example, a measurement of 5 °C at Lake Constance in the middle of August should be detected as an outlier, as it is 10 °C to 20 °C below the expected temperature for August. However, since temperatures of 5 °C are usual in January, the measurement is not conspicuous when performing straightforward outlier detection. For robust outlier detection on a temperature time series, the known temperature cycles such as the yearly temperature cycle (YTC) must be removed first. In our case, this is done by computing the mean annual curve from the raw time series and then subtracting it from the respective values of the time-series. The YTC-corrected time series can then be analyzed. For the outlier detection, we have used a robust outlier detection method described in Leys et al. [35] using the mean absolute deviation and the median.

For the analysis of trends, the raw time series is not suitable, as neither the data density nor the overpass-time where the temperature is measured is stable through time. This can induce additional trends. In Figure 4, one can nicely see the evolving data density through the years: In the 1980s, only one satellite at the time was available, whereas in 2013 there were five. Also, the heterogeneity of overpass-times in general, as well as the shift in overpass time due to orbit drift of the individual satellites, can be also be seen in Figure 4.

To overcome the data density issue, monthly and seasonal averages of the YTC anomaly are computed. The YTC anomaly is computed with reference to lake-specific YTC issued form the whole period of available LSWT data (1981–2016). The averaging over months or seasons leads to a homogeneous distribution of the temperature data suited for trend analysis. By using the YTC anomaly instead of the LSWT, we aim at removing an eventual bias induced by uneven distribution of the measurements within the averaging period (month or season).

Finally, to solve the overpass-time problem, a diurnal temperature cycle (DTC) model was used to adjust each measurement to the daily mean temperature, as described in the following paragraphs.

3.3. Diurnal Temperature Cycle (DTC) Correction

The NOAA satellites are subject to relatively strong orbit drifts [36,37]. This means that through the years of operation, the overpass time at the same location is shifted. Figure 4 shows all scenes and recording times used during this study. One can see that the overpass time of each satellite changes over the years of operation, leading to shifts of over five hours. As lake surface water temperature is changing throughout a day, an overpass time at 1 p.m. would lead to higher average temperatures than an overpass time at 8 am. Therefore, drifting orbit leads to an artificial trend in long-time records when looking at only one platform. When it comes to time series over several platforms, as is the case here, the impact is hard to predict. Some authors claimed that the impact is negligible (e.g., [32]); others correct for the DTC (e.g., [38]).

For this study, we used a DTC correction based on the DTC model described in Equations (2)–(5) after Duan et al. [39]. The model was first developed by Inamdar et al. [40] for land surface temperature (LST) and tested against in-situ data with various land cover types, including water (Lake Tahoe). Duan et al. [39] improved the model by deducing the width ω (Equation (4)) over the half-period of the cosine term from the thermal diffusion equation. This model was successfully implemented and used for daytime LSWT retrieval by Pareeth et al. (e.g., [38]).

$$\begin{gathered} T_s^d\left(t\right)=T_0+T_{a}cos\left(\frac{\pi }{\omega }\left(t-t_m\right)\right), t<t_s \\ {T}_s^n\left(t\right)=T_0+\delta T+\left[T_{a}cos\left(u\right)-\delta T\right]\frac{k}{\left(k+t-t_s\right)}, t\ge t_s \end{gathered}$$

(2)

with

$$k=\frac{\omega }{\pi }\left[tan{\left(u\right)}^{-1}-\frac{\delta T}{T_a} sin{\left(u\right)}^{-1}\right]$$

(3)

$$\mathsf{\omega }=\frac{4}{3}\left(t_m-t_{sr}\right)$$

(4)

$$u=\frac{\pi }{\omega }\left(t_s-t_m\right)$$

(5)

$T_s^d$ and $T_s^n$ are the day and night time surface temperatures, $t$ is the time, $T_0$ the residual temperature at sunrise, $T_a$ is the temperature amplitude, $\delta T$ is the difference in temperature between $T_0$ and $T_{\left(\mathrm{t}\to \infty \right)}$, $\mathsf{\omega }$ is the width over the half-period of the cosine term, $t_{sr}$ is the sunrise time, $t_m$ is the time where daily maximum temperature is reached, $t_s$ is the time where free attenuation starts, and $k$ is the attenuation constant.

During the here-presented study, instead of adjusting the satellite derived LSWT to a temperature equivalent of a standard acquisition time (12:00 UTC), we adjusted the temperature to the daily average temperature, computed from the modeled DTC curves. For each lake and for each month, typical DTC anomaly curves were computed by fitting the model to the available satellite LSWT data. The final LSWT was then adjusted to the daily average temperature by subtracting the temperature anomaly. The use of this approach improved the validation results within the here-presented study (see Section 4.2).

3.4. Trend Analysis

The linear trends were computed using the homogenized time series at the 26 locations (Table 1) described in the previous chapter. Non-parametric trend analysis was performed using the Thiel–Sen slope to quantify the trend and the Mann–Kendal test to determine its significance.

Several studies do compute only summer trends (e.g., [17]). This has the advantage of not having to deal with ice onset and breakup cycle. However, Winslow et al. [41] showed that seasonal trends such as a summer trend do not well represent the effects of climate change on the lakes. During this study, trends for all four seasons, as well as annual trends, were computed. When it comes to the definition of the seasons, there are mainly two ways described in literature: (1) summer defined as June, July, August (JJA), e.g., [42,43,44]; (2) summer defined as July, August, September (JAS), e.g., [17,32,41,45]. However, in the case of the lakes studied here, the influence of spring warming on the JJA period and the influence of fall cooling on the JAS period are disproportionately high (see Figure 5). We concluded that both JAS and JJA were not suitable for “summer only trend” in a restricted sense. For this study, we overcame this issue by taking the middle way; thus, the periods start and end in the middle of the months: Summer is defined as the period from 15 June to 14 September (JJAS), autumn as from 15 September to 14 December (SOND), winter from 15 December to 14 March (DJFM), and, finally, spring from 15 March to 14 June (MAMJ). Choosing these periods, the characteristics of the seasons are more pronounced (spring = warming, fall = cooling, summer = hot, and winter = cold) looking at the annual temperature cycle (Figure 5).

4. Results

4.1. Validation of Retrieved Lake Surface Water Temperature

The main performance indicator for the retrieval validation during this study is the standard deviation (SD), whereas the root mean square error (RMSE) and the bias are considered as secondary performance indicators. As we are primarily interested in temperature patterns such as temporal trends or spatial gradients, the uncertainty between two distinct measurements in time or space (SD) needs to be as small as possible, whereas the systematic part of the uncertainty (bias) is of lower importance. In fact, the bias affects all measurements in the same way, and therefore does not mask patterns such as spatial gradients and temporal trends. For the LSWT retrieval, the variable part of the uncertainty (SD) is mainly induced by variability in atmospheric composition, weather conditions, small-scale lake dynamic processes, and sensor noise, whereas the systematic part of the uncertainty (bias) is induced by factors such as skin-bulk temperature difference, in-situ sensor location, or sensor-related issues.

Further, it must be mentioned that during validation, we are comparing two measurements of a physical quantity against each other. This means that none of the two measurements (in-situ and satellite) do represent the truth to the full extent. This got even worse when considering the fact that they technically do not measure exactly the same quantity: The in-situ temperature measurement deduces the water temperature at a specific point in space (location and depth), whereas the satellite measurement estimates the temperature of the uppermost skin layer, averaged over a certain spatial extent, which is defined by its resolution. The exact amount of uncertainty issued from the described difference in measurement technique and measured quantity could not be determined during this study. Also, the correction methods for the skin to bulk difference as proposed by Donlon [46] or Minnet [47] for SST retrieval did not improve the validation performance: The method, applied to LSWT, only affects the systematic bias and has no effect on the variance of the error.

In Figure 6, one can see the difference when validating the same satellite data with data from two different in-situ locations close to each other. This relativizes the significance of such validation results, even though they still give a general indication about how well the retrieval works. However, it reminds us to be careful with the in-situ data, and to not take it as the finite truth.

Table 4. Validation results for each in-situ station. The data is ordered by latitude from south to north. The measurement type describes the averaging of the in-situ measurement: daily averaged, hourly averaged, and instant measurements. Instant measurements are single measurements that were manually taken during measurement campaigns. As the exact timestamps for the instant measurements where not known at Lake Balaton (Kesztehly old) and Lake Trichonida, 12:00 local time is assumed.

Location	Type	Location		n Match	R	RMSE	BIAS ± SD
Lake Trichonida	Instant *	38.57°N	21.56°E	16	0.987	1.59	−1.23 ± 1.04
Geneva INRA	daily	46.38°N	6.47°E	2453	0.985	1.08	0.26 ± 1.05
Geneva Buchillon	hourly	46.46°N	6.40°E	2373	0.981	1.29	−0.52 ± 1.18
Balaton Kesztehly old	Instant *	46.74°N	17.29°E	114	0.982	1.54	−0.57 ± 1.44
Balaton Keszthely new	hourly	46.76°N	17.27°E	1050	0.966	1.45	0.69 ± 1.27
Balaton Siofok	daily	46.90°N	18.05°E	2834	0.985	1.47	−0.54 ± 1.36
Constance Bregenz new	hourly	47.51°N	9.74°E	1723	0.990	1.01	0.36 ± 0.95
Constance Bregenz old	daily	47.51°N	9.74°E	719	0.981	1.37	0.50 ± 1.28
Constance Lindau	hourly	47.54°N	9.68°E	1755	0.976	1.53	0.63 ± 1.40
Müritz	instant	53.42°N	12.70°E	66	0.989	0.82	−0.07 ± 0.82
Peipsi Mehikoorma	daily	58.24°N	27.48°E	2102	0.984	1.25	0.44 ± 1.17
Pyhäjärvi Säkylä	daily	61.11°N	22.17°E	1749	0.975	1.26	0.36 ± 1.21
Päijänne Päjätsalo	daily	61.29°N	25.57°E	1017	0.947	2.06	1.28 ± 1.61
Päijänne Kalkkinen	daily	61.47°N	25.59°E	722	0.964	1.49	−0.61 ± 1.36
Pyhäselka Joensuu	daily	62.60°N	29.77°E	1893	0.962	1.46	0.11 ± 1.46
Lappajärvi	daily	63.26°N	23.63°E	2044	0.949	1.57	0.25 ± 1.55
Pielinen	daily	63.54°N	29.13°E	1707	0.938	1.83	0.45 ± 1.77
Oulujärvi	daily	64.40°N	27.17°E	1973	0.926	1.93	0.56 ± 1.84
Inarijärvi	daily	68.90°N	28.10°E	1207	0.896	2.60	1.65 ± 2.01

* No exact timestamp available (day only). 12 a.m. local time is assumed.

shows the validation results at all in-situ stations. There seem to be two main factors influencing the validation performance. First, there is the latitude, as it seems that with higher cloud cover rate and longer periods of ice cover during the cold season, also the quality of the retrieval decreases. On the other hand, the data type of the in-situ data seems to play an important role too: There are three types of in-situ data: daily averaged, hourly averaged, and manually taken instant measurements. The validation against instant and hourly averaged measurements yields the best results, as it allows for close temporal matchups between satellite overpass and in-situ measurements (matchup window ±1 h).

Amongst the instant measurements, the validation at Müritz yields the best results (SD ± 0.82 K, Figure 7a), but due to its sporadic sampling frequency only 66 matchups where available. The performance of the instant measurements at Lake Balaton (SD ± 1.44 K, Figure 7b) and Lake Trichonidas (SD ± 1.04 K) is lower, but at these locations the exact in-situ measurement timestamp was unknown. The validation with hourly averaged measurements has considerably higher number of matchups, and the performances reach from ±0.95 K (Lake Constance, Bregenz new, Figure 6a) to ±1.40 K (Lake Constance, Lindau, Figure 6b). The validation with daily averaged in-situ data shows, generally, a lower performance than with instant and hourly measurements. The SD reaches from very low ±1.07 K (Geneva, INRA Figure 8a) to rather increased values of ±2.05 K at Lake Inarijärvi (Figure 8b).

4.2. Validation of the DTC Correction and Trends

For the validation with daily in-situ data, the available satellite LSWT is compiled to daily average temperatures. To reduce the impact of the daily temperature cycle (DTC), a DTC correction is applied to the satellite retrievals. The validation with DTC-corrected data results in a decrease of SD by 0% to 15.6% depending on the lakes, and in an average decrease of 6% or 0.11 K, respectively (

Table 5. Impact of the DTC correction on the validation with daily in-situ data. Whereas for Lake Balaton Siofok, the improvement is neglectable, a considerable performance improvement is observed, for example, at Lake Constance, Bregenz old.

Location	BIAS ± SD		SD Change
	Raw LSWT	DTC Corr.
Balaton Siofok	−1.05 ± 1.37	−0.54 ± 1.36	−0.1%
Päijänne Päjätsalo	1.11 ± 1.63	1.28 ± 1.61	−1.4%
Peipsi Mehikoorma	0.01 ± 1.22	0.44 ± 1.17	−3.5%
Pyhäjärvi Säkylä	0.12 ± 1.26	0.36 ± 1.21	−4.2%
Pyhäselka Joensuu	−0.07 ± 1.52	0.11 ± 1.46	−4.3%
Päijänne Kalkkinen	−0.75 ± 1.42	−0.61 ± 1.36	−4.7%
Pielinen	0.31 ± 1.87	0.45 ± 1.77	−5.1%
Geneva INRA	−0.09 ± 1.11	0.26 ± 1.05	−5.5%
Lappajärvi	0.13 ± 1.67	0.25 ± 1.55	−7.5%
Oulujärvi	0.52 ± 2.00	0.56 ± 1.84	−7.9%
Inarijärvi	1.59 ± 2.37	1.65 ± 2.01	−15.5%
Constance, Bregenz old	−0.01 ± 1.51	0.50 ± 1.28	−15.6%

).

The trends from the in-situ stations with a period covered of at least 25 years were compared to the trends computed from LSWT retrievals. This was done for both the DTC-corrected satellite data and for the non-adjusted satellite data.

Table 6. Comparison of trends from satellite-derived LSWT with trends from in-situ data records. Only in-situ stations with data covering at least 25 years were used. The difference between DTC-adjusted trends and in-situ trends, as well as the difference between raw LSWT trends and in-situ trends, is shown in the last two columns.

Locations	Yr	In-Situ Trend		DTC Trend		RAW Trend		Diff. with In-Situ
		senn	p	senn	p	senn	p	DTC	RAW
Inarijärvi	31	0.42	0.000	0.45	0.000	0.58	0.000	0.03	0.16
Pyhäselka Joensuu	32	0.58	0.000	0.59	0.000	0.56	0.000	0.00	0.02
Peipsi Mehikoorma	26	0.44	0.000	0.44	0.003	0.45	0.004	0.00	0.01
Lappajärvi	31	0.45	0.000	0.50	0.000	0.22	1.000	0.05	0.23
Pielinen	31	0.50	0.000	0.64	0.000	0.53	0.000	0.14	0.03
Oulujärvi	31	0.43	0.000	0.49	0.000	0.59	0.000	0.06	0.15
								0.05	0.10

shows results from this comparison. At three of the six locations, the DTC adjusted trends are considerably closer to the in-situ trends (Inarijärvi, Lappajärvi, and Oulujärvi). At two locations (Pyhäselka and Peipsi), the DTC adjustment has almost no effect, and the LSWT trends remain close to the in-situ trends. At lake Pielinen, however, the DTC adjustment has a negative effect: the raw trend is in fact closer to the in-situ trend.

4.3. Trend Analysis

The yearly and seasonal trends at all studied locations are shown in

Table 7. Yearly and seasonal trends at all studied locations ordered from South to North. The covered period is slightly below the maximum of 36 years, depending on availability of cloud-free scenes at each location. The grayed-out values are trends (K/dec) in which Kendall p-value > 0.They are not retained to compute the average trends showed on the last line. The strongest seasonal trend of each lake is marked with an asterisk. Blank values signify that there were not enough values to compute trends, mostly due to ice cover during the cold period.

Location	Period	Yearly		Winter		Spring		Summer		Autumn
	Yrs	senn	p	senn	p	senn	p	senn	p	senn	p
Djorf Torba	34.9	0.17	0.03	−0.08	0.72	0.07	0.67	0.47 *	0.00	0.20	0.23
Lake Trichonida	35.5	0.21	0.00	0.40 *	0.00	0.18	0.02	0.18	0.00	0.16	0.16
Lake Ohrid	35.5	0.29	0.00	0.44 *	0.00	0.18	0.10	0.32	0.00	0.23	0.03
Lake Bolsena	35.5	0.28	0.00	0.46 *	0.00	0.37	0.00	0.25	0.02	0.09	0.25
Sanguinet	35.4	0.21	0.00	0.11	0.43	0.41 *	0.00	0.22	0.03	0.11	0.37
Lake Yalpuh	35.4	0.36	0.00	0.15	0.55	0.43	0.05	0.63 *	0.00	0.23	0.09
Lake Garda	35.4	0.37	0.00	0.53 *	0.00	0.37	0.00	0.29	0.00	0.30	0.01
Lake Geneva	35.5	0.33	0.00	0.44 *	0.00	0.38	0.01	0.26	0.02	0.29	0.03
Lake Balaton	35.4	0.34	0.00	0.17	0.57	0.54 *	0.00	0.53	0.00	0.12	0.29
Lake Constance	35.5	0.39	0.00	0.50 *	0.00	0.36	0.00	0.32	0.00	0.35	0.00
Kiev Reservoir	35.4	0.30	0.00	0.01	0.66	0.24	0.21	0.45	0.00	0.46 *	0.03
Müritz	35.4	0.39	0.00	0.10	0.55	0.69 *	0.00	0.30	0.04	0.38	0.02
Lough Neagh	35.4	0.33	0.00	0.40	0.02	0.45	0.00	0.04	0.83	0.46 *	0.00
Vättern	35.5	0.49	0.00	0.46	0.00	0.57 *	0.00	0.32	0.01	0.56	0.00
Lake Peipus	35.3	0.24	0.01	0.54	0.10	0.13	0.53	0.25	0.08	0.32 *	0.05
Vänern	35.4	0.42	0.00	0.35	0.03	0.36	0.04	0.31	0.02	0.54 *	0.00
Pyhäjärvi	35.3	0.54	0.00	0.57	0.05	0.47	0.05	0.57 *	0.00	0.54	0.00
Lake Päijiänne	35.4	0.47	0.00	0.76	0.18	0.45	0.04	0.50 *	0.00	0.49	0.00
Lake Onega	35.2	0.76	0.00			0.84	0.00	0.85 *	0.00	0.63	0.00
Pyhäselkä	35.3	0.55	0.00			0.54	0.01	0.53	0.01	0.57 *	0.00
Lappajärvi	35.2	0.45	0.00			0.40	0.16	0.45	0.01	0.57 *	0.00
Pielinen	35.3	0.48	0.00			0.39	0.26	0.44	0.01	0.61 *	0.00
Thingvallavatn	35.3	0.38	0.00			0.49 *	0.00	0.35	0.01	0.49 *	0.00
Oulujärvi	35.3	0.53	0.00			0.52	0.03	0.45	0.00	0.68 *	0.00
Lokka Reservoir	35.3	0.51	0.00					0.52	0.00	0.89 *	0.00
Inarijärvi	35.3	0.50	0.00					0.46	0.01	0.82 *	0.00
AVG (p < 0.05)		0.40		0.44		0.47		0.41		0.51

* Season with the highest trend per lake.

. Only significant trends (Kendall p-value < 0.05) were retained for further analysis.

All annual trends where significant and positive, thus indicating a warming of the lakes water temperature (Table 7). The strongest annual warming trend is observed at Lake Onega with 0.76 K/dec, whereas the weakest trend is to be found at Djorf Torba Reservoir with 0.17 K/dec. The trend at Lake Onega is remarkable, as there is no similar trend within the studied lakes; the next highest trend is 0.54 K/dec for lake Pyhäjärvi. The trend at Djorf Torba has to be looked at with care, as there has been a strong siltation going on since the drought period in the 1990s: The volume loss was estimated at 35% in 2014, with the siltation still going on at a rate of 0.9% per year [48].

Figure 9 shows details of the temperature anomaly in annual resolution of lakes representing geographical extremes (N-S, W-E). The significance of a long time series including years from 1981 onwards is clearly shown with increased frequency of blue bars during the earlier years (1981–1995) in comparison to the latest years with higher temperatures. Especially the warm period of the last five years is obvious, independently of the geographical location of the lakes. However, these high temperatures are more pronounced in northern and eastern areas, whereas the lake in Iceland (Thingvallavatn) shows reduced temperature increase, probably due to the oceanic influence. Also, one of the southern lakes (Trichonida, Figure 9) shows a reduced increase of temperature, which is probably also influenced by the proximity to the Mediterranean Sea. In all studied lakes, we can identify spatial gradients of the magnitude of temperature trends: The trends are increasing from south (Lake Trichonida 0.21 K/decade, Figure 9a) to north (e.g., Inarijärvi 0.50 K/decade, Figure 9b), and also from west (e.g., Thingvallavatn 0.38 K/decade, Figure 9c) to east (e.g., Lake Onega 0.76 K/decade, Figure 9d). The gradient is strongest on the diagonal from the south-west to north-east, whereas the magnitude of the trends shows a significant correlation with the latitude (R² = 58%, Figure 10a) and a rather weak correlation with longitude (R² = 24%, Figure 10b).

When looking at the strongest seasonal trends in Table 7, one can see that at lower latitudes the winter trends and some spring trends are dominant. This is typically the case for the subalpine lakes such as Lake Geneva, Lake Garda, and Lake Constance, which all show dominant winter trends. When moving north, first the spring (e.g., Balaton, Müritz, and Vättern) and then the summer trends (e.g., Pyhäjärvi, Päijänne, and Onega) start to dominate. Finally, for the highest latitudes, the autumn trends are clearly dominating (all lakes north to Lake Onega).

Looking at the north–south gradient of the seasonal trends (Figure 11), we see that winter trends seem almost unaffected by the latitude in contrast to the other seasonal trends. Especially, autumn trends are very sensible to latitude: At low latitudes, temperature trends in autumn are less pronounced than winter or spring trends, whereas at northern latitudes they are the strongest seasonal trends.

The West–East gradient of the seasonal trends is less significant, and probably influenced by the transition from oceanic to continental climate. However, this could not be investigated in more detail with the available data.

5. Discussion

The validation shows that with careful processing, the difficulties related to the AVHRR sensor such as orbit drift and small georegistration errors and others could be overcome. The importance of the in-situ measurements has been emphasized during the validation, and it was shown that the quality of the in-situ measurement is significantly influencing the performance of the validation. This indicates that the satellite retrieval has reached an accuracy level at which more accurate in-situ measurements are needed to further assess the quality of the LSWT retrieval. Also, the derived trends could be well validated with the in-situ measurements.

Concerning the significance test of the trends, only significant trends (Kendal p-value < 0.05 Table 7) were retained for further analysis. In contrast to the yearly trends, not all seasonal trends fulfilled this criterion. In the case in which the Kendall p-value is above 0.05, there are two plausible interpretations: (1) There is a strong trend, but due to elevated variance in the data and/or due to short periods and low amount of data, the trend cannot be determined with confidence. With this, erroneous estimates of trends can be detected (e.g., winter trend at Lake Päijiänne Table 7); (2) There is no trend or a weak trend in the data, and the p-value is increased, because we are close to the null-hypothesis. In this case, we eliminate a potentially correct value of weak trend close to 0 K/decades such as winter and autumn trends from lake Sanguinet (Figure 12). Thus, it is important to emphasize that with the used method, we can detect significant trends, but we cannot detect trends near 0 K/decade with certitude. This explains the discrepancy when comparing the annual trend with the valid seasonal trends from Table 7, in which our intuition is that the average of the seasonal trends is more or less equal to the yearly trend.

Trends are sensitive to different factors such as the period in which a trend is determined, the homogenization of the data, and the type of trend (seasonal, yearly, etc.). Therefore, the comparison with trends from individual lakes in other studies is not meaningful. However, on the other hand, more general statements can be looked at in detail. For example, Schneider and Hook [17] mention elevated summer trends in Northern Europe of up to around 0.8 K/dec. Our results suggest the same tendency of stronger trends in north-east Europe, whereas the strongest trends are rather to be found in autumn than in summer. The approximatively 0.8 K/dec in summer is only reached in Lake Onega, whereas at Inarijärvi and at the Lokka Reservoir such elevated trends are occurring during autumn.

Another observation from Schneider and Hook [17], that lakes in Northern Europe are warming more rapidly than local air temperature, was confirmed (see e.g., Lake Vättern Figure 13a). According to our results, this phenomenon is not limited to northern Europe: it can also be observed at some lakes in central Europe such as Lake Constance (Figure 13b). However, the difference in trends between local air temperature and LSWT is reduced when moving south towards warmer lakes. Additionally, for Lake Bolsena (Figure 13c), the air temperature trend tends even to be more pronounced than the LSWT trends. However, the reasons for the strong LSWT trends compared to local air temperature were not further investigated here.

The seasonal trends we detected show a high intra-seasonal variability within the same lake (Table 7), and therefore we can state that none of the seasons are representative for the overall trend within a lake. This is perfectly aligned with the main statement of Winslow et al. [41], who claim that summer trends do not fully represent the effects of climate change on lake temperatures. When we compare the yearly temperature cycle deduced from the first and the last 15-year period of the dataset, we can observe that the temperature evolution is very heterogeneous over the year, and also very different from lake to lake. For example, at Lake Bolsena and Lake Constance (Figure 14a,b), we observe significant warming during winter and spring. However, then there are also some short periods in which the recent years were colder than the first years of the dataset (mid-summer at Lake Constance, end of summer and start of autumn at Lake Bolsena). In contrast to this, at Lake Vättern and Lake Onega (Figure 14c,d), a moderate to significant warming during the whole year is observed, whereas this is more pronounced in summer and autumn. Lake Onega shows a significantly earlier and more pronounced warming in spring and during summer, resulting in a remarkable summer temperature peak that is higher and earlier than in the first years of the dataset.

Further, also within the same water body, trends can change significantly. This result is in agreement with the results of [16], in which intra-lake heterogeneity of trends was studied in large lakes. At lake Peipus (Figure 15), the yearly trend in the open waters in the north indicates a warming of 0.24 K/decade, whereas the warming rate almost doubles (0.44 K/decade) within the narrower part of the lake around Mehikoorma, in which also an in-situ station is located. This shows that the location where the temperature is sampled has an impact on the resulting trends. However, the example at lake Peipus shows us also that trends computed with the data of an in-situ measurement taken at a specific location are not necessarily representative of the temperature trend of the entire lake.

6. Conclusions

During this study, we could show that the AVHRR LAC data (1 km) archive at University of Bern can be used to yield good quality LSWT of an exceptional long-time span of over 35 years. This offers the possibility to study lake dynamics related to climate change. The quality of the retrieved LSWT could be confirmed by validation considering in-situ measurements, but it could also be shown that the quality of the validation data is essential.

A homogeneous LSWT dataset from 26 European lakes covering the period from 1981 to 2016 was created. With this, we were able to compute LSWT trends and to compare average yearly temperature cycles of different periods, of different geographic locations.

For all locations, significant LSWT trends indicating a warming were found. A gradient of increasing warming trends was revealed when moving from the south-west to the north-east of Europe. Whereas winter trends are not significantly affected, autumn trends are very sensitive to the location of the lake, especially to the latitude. Looking at the dominant seasonal trends per lake, we observe a shift from stronger winter trends in the south to dominant autumn trends in northern lakes, passing by a mixture of strong spring and summer trends in-between.

Generally, the strong interannual heterogeneity of the trends was detected in all lakes, showing that using only a defined period of the year to determine the temperature trend is not representative. It could be confirmed that the different seasons contribute unequal to the yearly trends, and that trend analysis based on summer seasons only does not show the complete reality.

Also, intra-lake heterogeneity of trends was detected. The location of measurements influences the trend directly, implying that a one-to-one comparison of trends between in-situ and satellite data has to be handled with care. We could show that temperature trends computed with data from in-situ measurements from individual stations are not suited to represent trends of an entire lake. This emphasizes the importance of satellite-based lake monitoring, as it adds significantly to the understanding of the big picture.

The comparison with local air temperature trends showed that many lakes have a stronger warming rate than local air temperature. This effect was shown to be more significant in colder lakes and in higher latitudes, e.g., Lake Vättern shows a stronger increase of LSWT than air temperature, but for Lake Bolsena in Italy the LSWT trend is even slightly weaker than for local air temperature.

The dataset generated within this study will be further investigated, and it has the potential to reveal further findings about the effects of climate change on European lakes. An extension of the dataset to more lakes in Europe would also be of great use for further investigations, especially regarding the east–west gradient of temperature trends, for which evidence was found. Also, a closer look at the influence on temperature trends of morphologic parameters such as the average depth of lakes could be the topic of further investigations.

For an extension of the dataset to other lakes in Europe, the size of the lakes is the only limiting factor, whereas for an extension on a global scale, the availability of AVHRR LAC data would be limiting. The application of the method to more recent sensors with higher resolution would be interesting for lake dynamic studies, but it is limited to a more recent time-period due to the availability of adequate satellite data.

Using open street map data for the generation of the LWM (land water mask) worked well. However, in the case of a larger number of lakes, and when taking in account the changes in time of global surface water occurrence, a global water mapping dataset such as presented in Pekel et al. 2016 [49] could be a more robust alternative for the generation of the LWM.