Synergy between satellite altimetry and optical water quality data towards improved estimation of lakes ecological status

: European countries are obligated to monitor and estimate ecological status of lakes under European Union Water Framework Directive (2000/60/EC) for sustainable lakes’ ecosystems in the future. In large and shallow lakes, physical, chemical, and biological water quality parameters are inﬂuenced by the high natural variability of water level, exceeding anthropogenic variability, and causing large uncertainty to the assessment of ecological status. Correction of metric values used for the assessment of ecological status for the effect of natural water level ﬂuctuation reduces the signal-to-noise ratio in data and decreases the uncertainty of the status estimate. Here we have explored the potential to create synergy between optical and altimetry data for more accurate estimation of ecological status class of lakes. We have combined data from Sentinel-3 Synthetic Aperture Radar Altimeter and Cryosat-2 SAR Interferometric Radar Altimeter to derive water level estimations in order to apply corrections for chlorophyll a, phytoplankton biomass, and Secchi disc depth estimations from Sentinel-3 Ocean and Land Color Instrument data. Long-term in situ data was used to develop the methodology for the correction of water quality data for the effects of water level applicable on the satellite data. The study shows suitability and potential to combine optical and altimetry data to support in situ measurements and thereby support lake monitoring and management. Combination of two different types of satellite data from the continuous Copernicus program will advance the monitoring of lakes and improves the estimation of ecological status under European Union Water Framework Directive.


Introduction
Consistent monitoring and estimation of the ecological status of lakes is essential for having an overview of current conditions and prediction of water quality (WQ) [1,2]. As lakes regulate the climate and carbon cycle, provide habitat for biota, drinking water, and other ecosystem services, it is important to monitor their ecological status and indications of their decline. To implement water management for sustainable environment, it is essential to evaluate the influence of human activity to lake ecosystems and also take into account variability of natural factors [3,4]. The European Union (EU) has established Water Framework Directive (2000/60/EC) [5] (WFD) for all European countries to achieve at least "good" ecological status class in lakes through consistent monitoring and water policy by the latest 2027 [6,7]. In total, five ecological quality classes from "high" to "bad" have been determined for every defined natural water body type, according to different biological (e.g., phytoplankton), physico-chemical (e.g., transparency) and hydro-morphological (e.g., depth variations) quality elements [5]. The main emphasis is given to the biological quality elements while the physico-chemical and hydro-morphological quality elements only affect the decision indirectly through their effect on the state of biological elements [5].
Traditional in situ measurements (e.g., collection of water samples, laboratory analyses, gauging stations, etc.) tend to be money and time consuming, and tend to take a lot of resources [8]. It is difficult to collect data from remote areas, where infrastructure is missing and access to the water body is complicated [9]. In situ gauging stations provide more frequent data than traditional measurements, but transmission of the data to the further workflow is rather slow and time consuming [10][11][12]. These kinds of systems are usually expensive to buy and difficult to maintain, therefore the number of gauging stations has decreased in recent years [11][12][13][14]. As lakes are mostly categorized as heterogeneous bodies, the biggest disadvantage of traditional in situ measurements is providing only point measurements, which do not represent spatial and temporal variations of the entire water body [1,11,15]. The combination of Earth Observation (EO) and in situ measurements complement each other by fulfilling their shortages and also helps to optimize and improve the selection of sampling locations and surveying times [16][17][18]. Earth Observation gives opportunity to have large scale data for spatial analysis and synchronized view of group of lakes [8,19]. In addition, digital data forms accelerates data reading and processing time and helps to produce historical and long time series [15]. The European Space Agency's (ESA) Copernicus program [20] provides full, free, open access, high quality near-realtime EO data on a global level for all users. Since the launch of the first Sentinel mission satellite in 2014, the purpose of the program is to launch nearly 20 satellites by the year of 2030 to provide long-term time series for various applications, e.g., environmental protection, climate change monitoring, and civil security [21]. The continuity of data allows to develop and advance methods to better integrate EO data into environmental monitoring, management, and policy.
As the satellite era has lasted almost 40 years by now, many algorithms and methods have been developed for different domains in lake research [10,[22][23][24][25][26][27][28][29][30][31]. As the WFD requirements are obligatory to every EU country, it has been indicated that the in situ methods are not feasible to cover the monitoring expectations, therefore many EO-based applications have been developed for implementation of requirements of WFD [32]. According to the WFD, phytoplankton is the key parameter of estimating ecological status of water [33][34][35][36], which also indicates eutrophication status of the lake [37]. Estimation of the concentration of the photosynthetic pigment chlorophyll a (chl a) has been one of the main goals in all ocean color missions starting from the first National Aeronautics and Space Administration (NASA) Coastal Zone Color Scanner on board Nimbus-7 in 1978. Additionally, water transparency is one of the key indicators for lake's ecological status class according to the WFD [38][39][40], where also EO-based assessment methods have shown the suitability to obtain additional information [7,25,32,41]. Covering various ocean color missions, the approaches and algorithms have been developed to provide trophic status based estimations [31,[42][43][44][45][46][47][48]] to obtain chl a and transparency estimates. However, estimation is based only on certain parameters (e.g., chl a and transparency), which do not take into account any external factors [34,35,38,49]. Ranges of different ecological status classes of lakes may be very narrow, which is the main reason why accurate WQ parameter estimation is essential. Water level (WL) estimation, using satellite altimetry data over inland water bodies, has become more advanced during recent years [10,16,17,[50][51][52], which allows the monitoring of its influence on the physical environment, which in turn changes the biota and whole lake ecosystem [53], especially in the case of shallow lakes.
In general, the ecological status of lakes is affected by anthropogenic factors, such as agricultural pollution by nutrients and industrial waste and various chemical elements [54], nevertheless natural factors, such as diurnal and seasonal changes, and prolonged changes in atmospheric circulations should be taken into account as well [55][56][57]. Water level is one of the natural factors, which affects WQ directly, worsening WQ parameters in low WL periods [58,59]. Especially in large, shallow lakes, the influence of lake depth variations on biota is remarkable, due to large seasonal and inter-annual WL fluctuations. The periods of extremely low WL are the most critical for the ecological status in such lakes, strengthening sediment resuspension, increasing concentrations of ions, phytoplankton biovolume and chl a concentration [60]. In contrast, the phytoplankton biomass (FBM) is significantly lower due to light limitation in years of high WL [61]. Such high variation in biological quality indicators, caused by natural factors, may over-mask human induced effects on a water body [62]. Based on the long-term in situ data, Tuvikene et. al. [62] showed the significant influence of WL while estimating the ecological status class according to the chl a, FBM and transparency in case of large shallow lakes. It was shown that, while phytoplankton species composition and traditional WQ parameters indicated different ecological status class, correction of metrics for WL enabled more consistent assessment of the ecological status of the lake.
As assessment of ecological status of lakes is currently parameter-based, thus accounting for the external natural factors is needed for more precise estimation. In the increasing constellation of ocean color satellites, carrying passive optical sensor usable for mapping various optical WQ parameters, and altimeters for estimating WL, the synergy of data could potentially improve the ecological status assessment over large shallow lakes. Therefore, the aim of this paper is to (1) create synergy between the optical and altimetry data in order to, (2) derive chl a, total FBM and transparency based ecological status classes taking into account the dependence of biological status indicators on the changing WL over large shallow lakes.

Study Sites Description
Data from two large and shallow lakes are used for this study. Võrtsjärv is situated in the central part of Estonia ( Figure 1) with surface area 270 km 2 (length 35 km, width~14 km). With the average depth of 2.8 m (the maximum depth is 6 m), the lake has strong WL fluctuation as the annual mean amplitude is 1.4 m and absolute range is 3.1 m, causing up to 3.2-fold difference in the lake volume [61]. Water level is generally characterized by a high period in spring, due to snow melting until end of April or beginning of May, a low period in winter and summer, and increasing WL from middle of the October due to autumn rainfalls. Võrtsjärv is a eutrophic lake with low water transparency (<1 m) [63]. According to the Estonian regulation of Water Act [64], Võrtsjärv belongs to independent ecological status class Type 6 ( Table 1) (non-stratified, light water color, low amount of chlorides, moderate hardness, surface area 100-300 km 2 ). Annually, national monitoring is performed at one station once in a month and from additional 2 or 9 stations in August [65].
Peipsi is a transboundary lake situated in the border of Estonia and Russia ( Figure 1) with surface area 3555 km 2 (length 152 km, width 48 km). The average depth is 7.1 m (the maximum depth is 15.3 m) and the annual mean amplitude of WL change is 1.5 m. Water level of Peipsi is characterized the highest in spring during the ice melting period and is the lowest in October [66]. It is a eutrophic lake, which consists of three parts: Peipsi sensu stricto (s.s., 2611 km 2 average depth 8.3 m), Lämmijärv (236 km 2 average depth 2.5 m) and Pihkva (708 km 2 , average depth 3.8 m). According to the Estonian regulation of Water Act [8], Peipsi belongs also to the independent ecological status class Type 7 (Table 1) (non-stratified, light water color, low amount of chloride, moderate hardness, surface area >1000 km 2 ) and is monitored 7 times per year (in March, when the lake is ice-covered, and from May to October) from 7 stations [65].
The both lakes are eutrophic and shallow, with similar WL fluctuation in terms of the annual mean (Võrtsjärv 1.4 m and Peipsi 1.5 m), but with different average depth (Võrtsjärv 2.8 m and Peipsi 7.1 m) leading to differences in the effect of WL. In addition, both lakes were considered unique, having no undisturbed reference sites for comparison, therefore, pressure-based approach, morphoedaphic index, and historical data were used to derive site-specific reference conditions and classification criteria for them [67]. (Võrtsjärv 2.8 m and Peipsi 7.1 m) leading to differences in the effect of WL. In addition, both lakes were considered unique, having no undisturbed reference sites for comparison, therefore, pressure-based approach, morphoedaphic index, and historical data were used to derive site-specific reference conditions and classification criteria for them [67].
Water level was measured by the Vaisala automatic weather station named as MAWS110, which sensor was located near the bottom of the lake and it measured the pressure of the water column, which is converted to the WL in centimeters [68]. Measurements were performed every hour, but WL was averaged over a day for the study. To derive WL above sea level (in Estonian Height system (EH2000), which corresponds to European Vertical Reference System), height of the sensor from sea level (national geoid named as EST-GEOID2017 (Estonian national geoid model)) [69] was added. Hydrometric stations were located near the shore or in the mouth of the river ( Figure 1). Water level data was available in Rannu-Jõesuu (hereafter Rannu) hydrometric station from 1927, in Mustvee and Praaga from 1921 and in Mehikoorma from 1949 until 2019.

Chl a, FBM, and SDD
Chl a and FBM samples of integral water at discrete depths (0.5 m, 1 m, etc., until 0.5 m from the bottom) were collected from different stations and analyzed in a laboratory from the 1960s. Chl a samples were filtered through Whatman glass-microfiber filters (GF/F) ø 25 mm filter and pigment extraction was done with 5 mL 96% ethanol [70]. Chl a was measured spectrophotometrically with Hitachi U-3010 (430-750 nm) and the concentrations of chl a were derived according to Jeffrey and Humphrey [71]. For determination of FBM, up to 400 counting units were counted and measured with an inverted microscope based on Utermöhl technique [72]. SDD values were measured with white Secchi disc from the 1950s.

Satellite Data
Sentinel-3A/3B (S3A/S3B) and CryoSat-2 (C2) data from 2016 to 2019 covering the period from May to October were used. Buffers of 1 km (Võrtsjärv) and 2 km (Peipsi) were created to reduce the possible errors in the vicinity of land in case of optical and altimetry data ( Figure 1). Depending on the cloud cover, 300 m resolution data from S3 Ocean and Land Colour Instrument (OLCI) is available daily in the study location. Surface elevations from the Synthetic Aperture Radar Altimeter (SRAL) sensor onboard S3A and S3B are available every 27 days for a given track. Lake Peipsi is measured by both S3A and S3B, hence, data are available from 2016 and onwards. Lake Võrtsjärv is only measured by S3B, hence data are available from the launch in 2018. C2, carrying Synthetic Aperture Radar/Interferometric Radar Altimeter (SIRAL) sensor, is operating in a geodetic orbit with a repeat period of 369 days. S3 and C2 tracks during the study period are shown on the Figure 1 in blue and gray, respectively.

Altimetry Data
To construct water surface elevations from S3A and S3B the ESA "Non Time Critical" Level-2 product "Enhanced measurements" was used. This product consists of 20Hz alongtrack measurements of, e.g., tracker range, satellite altitude, atmospheric and geophysical corrections, and waveforms. To identify the nadir water surface, the waveforms were retracked with a narrow primary peak threshold retracker [76]. For C2 we applied the ESA L1B product from baseline D. This product also contains 20Hz measurements of, e.g., waveforms, tracker range, satellite altitude and corrections. The same retracker was applied to the C2 waveforms.
The water surface elevations are derived from the following relation (Equation (9)) where h is the altitude of the satellite, and N is the geoid height, which here is the EST-GEOID2017 national geoid model. R is the range, i.e., the distance from the satellite to the surface which can be expressed as (Equation (10)) Here R tracker is the distance to the presumed surface measured by the onboard tracker and R retrack is the retracking correction. R atm represents the atmospheric corrections consisting of the wet and dry troposphere and ionosphere. R geo represents the geophysical corrections including solid earth tide, pole tide, and ocean loading tide. All corrections are here taken from the altimetry products.
The altimetry-based WL time series based for Peipsi and Võrtsjärv were reconstructed using the "R" package "tsHydro" available from https://github.com/cavios/tshydro (accessed on 14 October 2020). The software package is based on an implementation of a simple state-space model where the observations are assumed to be described by a random walk plus noise. The observation noise is described by a mixture of a Normal and the heavy tailed Cauchy distributions making it more robust against erroneous observation. In the solution, bias between the different missions is also estimated. The model is described in detail in [10]. Here, one time series based on all available altimetry data is constructed for each lake.

Data Correction for the Effects of Water Level
EO-derived WQ variables (chl a, FBM, SDD) were corrected for the effects of changing WL by the methodology developed by Tuvikene et al. [62]. The methodology is based on the long-term in situ data, where statistically significant relationships (p-value < 0.05) between the specific WQ variables and WL were found (Table A1). Corrected WQ variables were calculated using regression equations according to the long-term monthly mean WL and daily measured WL. If the WL of the sampling day corresponded to the long-term average for that period, no correction was needed, but if it differed, the correction either up or down was proportional to the deviation of the WL from the mean.

Inter-Annual Dynamics in Water Level
Comparison of long-term WL measured in three hydrometric stations (Mustvee, Mehikoorma, Praaga) in Peipsi ( Figure 2) show similar dynamics and very small spatial differences, e.g., the mean difference between Praaga and Mustvee is 0.02 m, and between Mehikoorma and Mustvee is 0.04 m. In both lakes, the higher and lower WL periods vary similarly and long-term monthly mean WL is highest in May and decreases towards October. While the annual difference can exceed the long-term monthly mean WL (Figure 2) in some years more than 1 meter, its effect is more pronounced in shallower Võrtsjärv (average depth 2.8 m) compared to Peipsi (average depth 7.1 m).
Water levels based on altimetry data and WL measured at the Mustvee hydrometric station show very similar dynamics ( Figure 3) with a bias of 0.4 m. Bias corrected Root Mean Standard Error (RMSE) is 0.07 m, whereas the Pearson correlation between WL measured in Mustvee hydrometric station, and derived from altimetry data, is 0.99. Comparison between altimetry data and Rannu hydrometric station in Võrtsjärv shows similar results, with a bias 0.4 m. Bias corrected RMSE is 0.08, whereas the Pearson correlation between WL measured in Rannu hydrometric station and WL, derived from altimetry data, is 0.98. The combination of altimetry data from C2 and S3 missions allows to obtain much more frequent time series, especially since 2018 due to the launch of S3B (Figure 3). October. While the annual difference can exceed the long-term monthly mean WL ( Figure  2) in some years more than 1 meter, its effect is more pronounced in shallower Võrtsjärv (average depth 2.8 m) compared to Peipsi (average depth 7.1 m). Water levels based on altimetry data and WL measured at the Mustvee hydrometric station show very similar dynamics ( Figure 3) with a bias of 0.4 m. Bias corrected Root Mean Standard Error (RMSE) is 0.07 m, whereas the Pearson correlation between WL measured in Mustvee hydrometric station, and derived from altimetry data, is 0.99. Comparison between altimetry data and Rannu hydrometric station in Võrtsjärv shows similar results, with a bias 0.4 m. Bias corrected RMSE is 0.08, whereas the Pearson correlation between WL measured in Rannu hydrometric station and WL, derived from altimetry data, is 0.98. The combination of altimetry data from C2 and S3 missions allows to obtain much more frequent time series, especially since 2018 due to the launch of S3B (Figure 3).  (5) to October (10). Long-term monthly mean WL is showed in solid line.

Inter-Annual Dynamics of Water Quality Parameters
Long-term in situ time series of chl a, FBM, and SDD show variation in seasonal dynamics and pronounced inter-annual differences in Võrtsjärv ( Figure 4A). For example, chl a shows increasing trend from 1995 until 2008, then it started to decrease again ( Figure 4A). FBM shows increasing trend from the early 2000 s. SDD shows seasonal variation from less than 0.5 to 1.5 m. As WQ parameters are in situ measured once during the vegetation period in Võrtsjärv station named as V2 (Figure 4B), the seasonal dynamics are better obtained with more frequent satellite data resulting on average 50 points, which doubled with S3B launch. This allows to acquire more representative data to analyze the WQ parameters over specific time period and have better knowledge of the seasonal and inter-annual differences.
The dynamics of satellite data corresponded well to in situ measured WQ parameters, showing similar dynamics during vegetation period with a bias for chl a 4.9 ± 11.8 mg/m 3 , for FBM 2.7 ± 14.2 g/m 3 , and for SDD was 0.3 ± 0.8 m. Chl a derived from OLCI corresponds to in situ data with Pearson correlation 0.7, whereas FBM 0.6 and SDD 0.4.
Depending on the monitoring system, some areas are not covered with sufficient monitoring frequency ( Figure A1

Inter-Annual Dynamics of Water Quality Parameters
Long-term in situ time series of chl a, FBM, and SDD show variation in seasonal dynamics and pronounced inter-annual differences in Võrtsjärv ( Figure 4A). For example, chl a shows increasing trend from 1995 until 2008, then it started to decrease again ( Figure  4A). FBM shows increasing trend from the early 2000′s. SDD shows seasonal variation from less than 0.5 to 1.5 m. As WQ parameters are in situ measured once during the vegetation period in Võrtsjärv station named as V2 ( Figure 4B), the seasonal dynamics are better obtained with more frequent satellite data resulting on average 50 points, which doubled with S3B launch. This allows to acquire more representative data to analyze the WQ parameters over specific time period and have better knowledge of the seasonal and inter-annual differences.  The dynamics of satellite data corresponded well to in situ measured WQ parameters, showing similar dynamics during vegetation period with a bias for chl a 4.9 ± 11.8 mg/m 3 , for FBM 2.7 ± 14.2 g/m 3 , and for SDD was 0.3 ± 0.8 m. Chl a derived from OLCI corresponds to in situ data with Pearson correlation 0.7, whereas FBM 0.6 and SDD  (Table A1). In Võrtsjärv, chl a was significantly related to WL (p < 0.05) in July, August, and October, FBM was significantly related in July, whereas SDD was related in all months (from May to October). In Peipsi s.s., WL was significantly related only with FBM in June and September and in Pihkva, FBM was related in August. In Lämmijärv, chl a was significantly related in July and August, FBM was related in June and SDD in June, July, and August. When the relationship was significant between WL and specific WQ parameter, correction was applied on a daily basis according to the WL value in respect to long-term monthly mean WL. When WL was higher than long-term monthly mean WL, WQ parameter values increased after correction, if the WL was lower, then WQ parameters were corrected to lower values. This changed the mean value of the month, as well as minimum and maximum values.

Relationships between WL and WQ parameters vary in both lakes
In Võrtsjärv, WL was higher than long-term monthly mean WL in 2016 (up to 0.3 m) and lower in 2017 until 2019, especially low in 2018 and 2019 towards autumn, when WL was lower by more than 0.5 m ( Figure 5). As WL was high in 2016, the mean value of chl a and SDD increased after applying the WL correction, whereas, in 2017 until 2019, mean values decreased. Due to the decrease in mean values, ecological status class of WQ parameters also changed. In 2017 and 2018, the ecological status class of chl a changed from "moderate" to "good", and in 2019 from "bad" to "moderate", whereas the ecological status class remained the same in 2016, although mean value of chl a increased. The ecological status class of SDD, changed from "good" to "very good" in 2018 and 2019, when the WL was the lowest compared to long-term monthly mean WL ( Figure 5).
In Peipsi, there was less significant relationships between WL and WQ parameters than in Võrtsjärv (Table A1). In Peipsi s.s., there was no statistically significant relationships between WL and WQ parameters, except for FBM in June and September, which only leads to a small change of mean values of FBM ( Figure A2 Peipsi s.s. FBM). As WL was lower than long-term monthly mean WL in 2016, 2018, and 2019, mean FBM value slightly decreased without changing the ecological status class. More significant relationships between WL and WQ parameters appeared in Lämmijärv ( Figure A2 Lämmijärv), which also changed the ecological status class of SDD. As WL was more than 0.4 m lower than long-term monthly mean WL in 2019, then ecological status class of SDD changed from "moderate" to "good", and the mean value of SDD in 2016 also decreased. In Pihkva, the amount of in situ data was mostly insufficient for statistical analysis (significant relationship was found only in August for FBM) and therefore, there were no changes in ecological status class estimations of WQ parameters ( Figure A2 Pihkva).
As Võrtsjärv and Lämmijärv are shallower than Peipsi s.s., statistically significant relationships between WL and specific WQ parameters tend to occur more frequently there on a monthly basis. More pronounced changes in WL (difference approximately from 0.4 m) from long-term monthly mean WL were influencing not only the mean values of WQ parameters, but even changed the ecological status class from one to another. Smaller changes in WL in respective to long-term monthly mean influenced only the mean values of WQ, not the ecological status class. However, it is still important, considering long-term monitoring of the lake. relationships between WL and specific WQ parameters tend to occur more frequentl there on a monthly basis. More pronounced changes in WL (difference approximatel from 0.4 m) from long-term monthly mean WL were influencing not only the mean value of WQ parameters, but even changed the ecological status class from one to anothe Smaller changes in WL in respective to long-term monthly mean influenced only the mea values of WQ, not the ecological status class. However, it is still important, considerin long-term monitoring of the lake.

Discussion
Long-term time series of WQ parameters in lakes are important to obtain information about the seasonal dynamics, their past and present situation and form the basis for future predictions [77,78]. The variety of measured parameters allow the partition of the influences by anthropogenic, as well as natural factors. Moreover, long-term time series help to understand the necessity for current management improvements for the future [79]. Water level as a natural factor, has an impact on lakes ecosystem and on WQ parameters, which in turn have a high impact on the ecological status of lakes [80,81]. Changes in the lake ecosystem influence biota existence and survival complexity in the lake, which affect ecosystem structure and future perspective [82][83][84]. Water level plays an important role in the lake ecosystem, especially in large shallow lakes [85], where monthly WL differences can be highly different from long-term observations, therefore removing WL influence from WQ parameters is beneficial for fair status assessment [62]. This study shows the possibility to create synergy between optical and altimetry data for lakes research, especially in terms of WFD reporting needs. The removal of one natural aspect will improve the estimation of the ecological status class of a lake, which will allow to keep the focus on the anthropogenic changes in the lake ecosystem.
The study shows high variations between seasonal and long-term WL in large and shallow lakes. Low WL periods vary with higher WL periods in both Võrtsjärv and Peipsi, which are caused by variation in precipitation, hydrological inflows, and evaporation in different years [86]. Long-term monthly mean WL decrease, from May towards October, characterizing higher WL after snow and ice melting, and lower WL in summer due to less inflow and higher rate of evaporation during the warm period [86,87]. Not only seasonal temperature and precipitation variability influences WL, but climate change has also direct effect to lake ecosystem [88][89][90][91], which has been even pointed out by The Intergovernmental Panel on Climate Change [92] referring to a global problem. The decrease in WL in some areas is often not compensated by the precipitation, which is caused by the higher annual mean air temperature, which, in turn, causes intensive evaporation and higher water temperature, causing a decrease in ice cover thickness and duration [93,94]. On the other hand, extreme events (flooding, rainfalls, storms, etc.) due to climate change are causing an abrupt increase in the WL in some areas, causing also more apparent WL fluctuations [81,86].
As WL is sensitive to external factors, it is important to monitor it continuously [86,90,[95][96][97]. As the number of in situ gauging stations is decreasing due to their cost and maintenance needs, the use of altimetry data becomes an important supplement to measure lake level variations. The use of altimetry has already provided promising results over lakes [10,11,16,50,[98][99][100][101][102][103][104][105]. This is confirmed in this study, where a good agreement between the altimetry-based WL and WL measured at hydrometric stations with a correlation up to 0.99 in both lakes, with RMSE for Võrtsjärv 0.07 and for Peipsi 0.08, respectively. By applying measurements from different altimetry satellites the amount of data is increased, hence, providing a better understanding of WL dynamics and its changes [106][107][108]. However, combining data located at different tracks may introduce errors due to unresolved geoid signals, or other effects that causes the lake surface so be spatially different such as wind or land signals near the shoreline, but could be combined when the spatial differences in WL are small inside the water body. For example, comparing WL in Peipsi hydrometric stations, seasonal dynamics of the WL were similar and differences remained between 0.02 and 0.04 m, which shows possible WL homogeneity in the water body. Here, S3 SRAL and C2 SIRAL data were combined to increase the data frequency, but other missions could also have been used. Additionally, S3 and C2 provide SAR mode data in the study area, which give higher along-track resolution than conventional altimetry, with dense C2 spatial tracks suitable for lakes. The continuation of S3 missions, and the recently launched S6 mission, ensure data for future monitoring, which allows to continue developing EO-based applications to monitor WL and its implication for ecosystems and their services.
The possibility to derive certain WQ parameters from satellites helps to collect data more frequently and from places, which are not easily accessible. Satellite data corresponded to in situ data showing correlation 0.7 for chl a and lower than 0.6 for FBM and SDD, whereas bias for chl a was 4.9 ± 11.8 mg/m 3 , for FBM 2.7 ± 14.2 g/m 3 , and for SDD 0.3 ± 0.8 m. More data from satellites gives an extra information to in situ data, which enables to detect seasonal dynamics and rapid changes in different parts of a lake, especially in terms of estimating chl a, FBM, SDD based ecological status class. In some cases, satellite data overestimated ( Figure 4 Võrtsjärv B SDD, Figure A1 Lämmijärv B SDD) or underestimated ( Figure A1 Peipsi s.s. B chl a) in comparison of in situ data. This can be due to spatial differences (point measurement vs. 300 m × 300 m pixel), measurement technique (integral sample vs upper water column captured by satellite), but still, both in situ and EO data show similar trends during the vegetation period and complement each other. Different WQ parameters have been used for estimating ecological status of lakes from satellite data by many authors already for a while [7,25,33,36,41,109], which confirms necessity of this work for WFD.
Based on the long-term in situ data, the study shows significant relationships (p < 0.05) between WL and WQ more frequently in Võrtsjärv than in Peipsi. It is caused by the fact, that Võrtsjärv is shallower than Peipsi and high seasonal WL amplitude affects its biota in greater extent. In addition, the relationships between WL and WQ parameters were more often significant in the shallower parts of Peipsi (e.g., in Lämmjärv). Correcting WQ parameters with WL results in the change of mean annual WQ parameter values and even in changes of their ecological status class. The impact was more pronounced in Võrtsjärv, where WL was much lower than long-term monthly mean WL (approximately from 0.5 m), which did not only decrease the mean value of WQ parameters, but even changed the ecological status estimation to a higher class. As WL has a high impact to WQ parameters and the whole ecosystem, it is encouraging to see that EO data can be used to remove WL effect from data of WQ estimates in large shallow lakes, thus it could be one of the research topics deserving more attention due to its applicability to other large and shallow lakes. Moreover, the variety of past, present, and future satellite missions allows to apply the demonstrated methodology to extend the temporal and also spatial extent of the results. This could be done by combining satellite sensors for WQ, e.g., Environmental Satellite (ENVISAT) Medium Resolution Imaging Spectrometer (MERIS), Moderate Resolution Imaging Spectroradiometer (MODIS) to estimate WQ parameters [35,40,[110][111][112][113][114], with Topex Poseidon and Jason mission satellites to estimate WL data [51,104,105,115,116]. Additionally, several databases like Hydroweb [105,117], Database for Hydrological Time Series of Inland Waters (DAHITI) [98,118], Global Reservoirs and Lakes Monitor database (G-REALM) of United States Department of Agriculture (USDA) [119] provide water level data of lakes, whereas Globolakes [120,121] and Climate Change Initiative (CCI) [122] provide data of WQ parameters from different satellites. However, Sentinel missions will provide data at least for 20 years, which gives a possibility to gather continuous data from now.
As sustainable lake ecosystems are essential in the present and in the future, it is important to ensure frequent monitoring system and to develop EO-based approaches to improve monitoring and assessment of the ecological status class of lakes, leading, finally, to the better management of lakes. In situ monitoring is still essential to gather data to monitor the ecosystem, but also it is basis for developing and improving new tools and approaches, while satellite data would support it with better spatial and temporal resolution for certain parameters. Synergy between altimetry and optical data, gives possibility to estimate WL, as well as WQ parameters for more precise analysis. Results of this study show that, WL can have a strong impact on WQ parameters in large shallow lakes, where amplitude in WL change is high. This effect can be corrected with the proposed method to improve the estimation of ecological status class of such lakes.

Conclusions
EU WFD obligates European countries to monitor WQ of lakes and estimate their ecological status according to reference condition. The ecological status of lakes is influenced by anthropogenic, as well as by natural factors. To decrease the effect of natural factors, WQ parameters should be corrected for WL in case of shallow lakes with changing WL for more accurate identification of human impact and estimation of ecological status. Traditional monitoring methods give valuable long-term time series, however often with constrained spatial and temporal coverage, while satellites provide more frequent data. Many water quality parameters (e.g., chl a, FBM, and SDD) and WL can be derived using EO based data. The development and validation of EO algorithms allow to create synergy between different types of satellite data which enables to analyze and assess the lake ecosystems and provide complementary information to improve the monitoring and management decisions.
Long-term in situ data show inter-annual and seasonal dynamics of indicators characterizing lake ecosystems. In large and shallow lakes such as Võrtsjärv and Peipsi, amplitude of seasonal WL changes is high (more than 1 m), which influences WQ parameters, especially in shallower Võrtsjärv. The results support previous studies [10,99,101,104] in the fact that WL can be estimated with high accuracy (bias 0.4 m, R 0.99) over lakes with EO data. Correcting WQ parameters with the WL data derived from satellite, the mean values of the parameters and resulting ecological status classes changed in some cases. The highest changes occurred in Võrtsjärv during the years with lower WL compared to the long-term monthly mean. In this case the correction of the values of WQ parameters changed the resulting ecological status to better quality class. Therefore, taking into account the variability of WL and correcting its effect in case of shallow lakes is important to decrease the influence of the natural factors on the assessment of the lakes ecological status class in terms of EU WFD.  Acknowledgments: We are grateful to anonymous reviewers for their feedback and suggestions that helped to improve the manuscript. This study was financially supported by (1) the European Regional Development Fund within National Program for Addressing Socio-Economic Challenges through R&D (project named as RITA1 grant number RITA1/02-52-01); (2) European Union's Horizon 2020 research and innovation program under grant agreement no. 730066; (3) Estonian Research Council grants PSG10 and PRG709. Authors thank people from Centre for Limnology providing in situ water quality data and Estonian Weather Service providing in situ water level data.

Conflicts of Interest:
The authors declare no conflict of interest. Appendix B Table A1. Relationship between in situ measured WL (cm) and water quality parameters from Appendix B Table A1. Relationship between in situ measured WL (cm) and water quality parameters from May (5) to October (10), where bold values indicate statistically significant relationships (p < 0.05), N indicates number of observations during the study period and "-" indicates not enough observations for statistical analysis.