Determining lake surface water temperatures worldwide using a tuned one-dimensional lake model (FLake, v1)

. A tuning method for FLake , a one-dimensional (1-D) freshwater lake model, is applied for the individual tuning of 244 globally distributed large lakes using observed lake surface water temperatures (LSWTs) derived from along-track scanning radiometers (ATSRs). The model, which was tuned using only three lake properties (lake depth, snow and ice albedo and light extinction coefﬁcient), substantially improves the measured mean differences in various features of the LSWT annual cycle, including the LSWTs of saline and high altitude lakes, when compared to the observed LSWTs. Lakes whose lake-mean LSWT persists below 1 ◦ C for part of the annual cycle are considered to be seasonally ice-covered. For trial seasonally ice-covered lakes (21 lakes), the daily mean and standard deviation (2 σ) of absolute differences between the

Abstract.A tuning method for FLake, a one-dimensional (1-D) freshwater lake model, is applied for the individual tuning of 244 globally distributed large lakes using observed lake surface water temperatures (LSWTs) derived from alongtrack scanning radiometers (ATSRs).The model, which was tuned using only three lake properties (lake depth, snow and ice albedo and light extinction coefficient), substantially improves the measured mean differences in various features of the LSWT annual cycle, including the LSWTs of saline and high altitude lakes, when compared to the observed LSWTs.Lakes whose lake-mean LSWT persists below 1 • C for part of the annual cycle are considered to be seasonally ice-covered.For trial seasonally ice-covered lakes (21 lakes), the daily mean and standard deviation (2σ ) of absolute differences between the modelled and observed LSWTs are reduced from 3.07 • C ± 2.25 • C to 0.84 • C ± 0.51 • C by tuning the model.For all other trial lakes (14 non-icecovered lakes), the improvement is from 3.55 • C ± 3.20 • C to 0.96 • C ± 0.63 • C. The post tuning results for the 35 trial lakes (21 seasonally ice-covered lakes and 14 non-icecovered lakes) are highly representative of the post-tuning results of the 244 lakes.
For the 21 seasonally ice-covered lakes, the modelled response of the summer LSWTs to changes in snow and ice albedo is found to be statistically related to lake depth and latitude, which together explain 0.50 (R 2 adj , p = 0.001) of the inter-lake variance in summer LSWTs.Lake depth alone explains 0.35 (p = 0.003) of the variance.Lake characteristic information (snow and ice albedo and light extinction coefficient) is not available for many lakes.The approach taken to tune the model, bypasses the need to acquire detailed lake characteristic values.Furthermore, the tuned values for lake depth, snow and ice albedo and light extinction coefficient for the 244 lakes provide some guidance on improving FLake LSWT modelling.

Introduction
The response of lake surface water temperatures (LSWTs) to climate is highly variable and is influenced by lake physical characteristics (Brown and Duguay, 2010).Some large lakes have been shown to alter the local climate.The extent of ice cover on lakes is considered to be a sensitive indicator of and also a factor in global climate change (Launiainen and Cheng, 1998).Changes in the length of the ice cover period affect local climatic feedbacks; for example, a shorter ice cover period allows for a longer time for surface heat exchange with the atmosphere (Ashton, 1986).This is of particular importance in areas where there is a high concentration of lakes, such as Canada (Pour et al., 2012).The fluxes of heat and moisture from the Great Lakes and the large Canadian lakes of Great Bear and Great Slave can impact the Published by Copernicus Publications on behalf of the European Geosciences Union.
mesoscale weather processes causing lake-effect storms, altering the local climate (Sousounis and Fritsch, 1994;Long et al., 2007).Shallow lakes, particularly those with a large surface area, such as Lake Balaton, are more sensitive to atmospheric events (Voros et al., 2010).
Reliable modelling of LSWTs can enrich our understanding of the highly variable dynamic nature of lakes.In this paper, a freshwater lake model, FLake (available at http: //www.flake.igb-berlin.de/sourcecodes.shtml), is tuned with along-track scanning radiometers (ATSR) Reprocessing for Climate: Lake Surface Water Temperature and Ice Cover (ARC-Lake) observations (MacCallum and Merchant, 2012) of 244 globally distributed lakes.FLake is a one-dimensional (1-D) thermodynamic lake model, capable of predicting the vertical temperature structure and mixing conditions of a lake (Mironov et al., 2010).The tuned model is expected to improve the LSWT representation of these lakes in FLake.
There have been some studies carried out that compare modelled LSWTs from FLake with LSWT observations on European lakes (Voros et al., 2010;Bernhardt et al., 2012;Pour et al., 2012).The findings from two of these three studies showed consistent differences between the modelled and observed LSWTs (overestimation of the open water LSWTs and underestimation of the ice cover period).Despite these differences, FLake is considered to be a reliable model for simulating LSWTs and lake ice phenology (Bernhardt et al., 2012).The overestimation of the open water LSWTs and underestimation of the ice cover period are consistent with preliminary trial work findings from this study, which included North American and European lakes.
It is the intention of this study to achieve a daily mean absolute difference (MAD) of ≤ 1 • C, between the tuned and observed LSWTs, across all lakes.For each lake, the daily MAD quantifies the mean absolute difference in modelled and observed daily LSWT value, averaged over the total number of days.The dispersion of the daily MADs are reported to 2 standard deviations (±2σ ) for each lake.When reported across a number of lakes, the daily MAD for each lake is averaged, with the dispersion of the MADs across individual lakes reported to ±2σ .
A daily MAD of ≤ 1 • C across all lakes, is possibly accurate enough for a global scale study.A lower daily MAD target may not be achievable as this study comprises of lakes with a wide range of geographical and physical characteristics.The effect of the tuning on the sub-surface temperature profile and on the depth of the mixed layer is not considered in this study.Many lake-specific properties can be considered in FLake.Preliminary model trial work was carried out on seven seasonally ice-covered lakes (deep and shallow), which had available lake characteristic data in the ILEC world lake database (ILEC, 1999) or LakeNet's Global Lake Database (LakeNet, 2003).Through this preliminary work, the lake-specific properties, which exerted the strongest effect on the modelled LSWTs, were selected.These properties are lake depth (d), snow and ice albedo (α) and light extinc-tion coefficient (κ).In the next part of the preliminary work, it was determined that the modelled LSWTs could be tuned to compare well with the observed LSWTs, by adjusting the values for these three properties: d, α and κ, herein referred to LSWT-regulating properties.On the basis of the preliminary findings, the trial work was performed on 35 lakes, prior to attempting to tune all 246 lakes.
An example of the preliminary trial work is shown in Fig. 1a, for Lake Athabasca, Canada (mean depth of 26 m).In this figure, a greater modelled α (higher reflectivity) results in a later ice-off date than the default model snow and ice albedo and is closely comparable to the observed ice-off date.In Fig. 1b, it is demonstrated that by using a shallower d than the mean depth of the lake, the ice-on day occurs earlier and corresponds more closely to the observed ice-on day.Lake depth is essentially being used as a means to adjust the heat capacity of the lake, exerting control over the lake cooling and therefore the ice-on date.The modelled LSWT is further improved by lowering the κ value (greater transparency).
The greater transmission of surface heat to the lower layer results in a lower (and more representative) maximum LSWT, Fig. 1b.The LSWTs modelled using the greater α, lower d and lower κ compare closely with the observed LSWTs, Fig. 1c.
In this study, the modelled mean differences for several features in the LSWT annual cycle are measured, quantifying the level of agreement with the observed ARC-Lake LSWTs, are the basis for selecting the tuned (optimal) LSWT-regulating properties (d, α and κ) for each lake.Lakes are divided into two distinct categories.Lakes with a lakemean LSWT climatology (determined using twice-a-month ARC-Lake full-year LSWT observations, 1992LSWT observations, /1996LSWT observations, -2011) ) remaining below 1 • C for part of the seasonal cycle are categorised as seasonally ice-covered lakes (160 lakes).All other lakes are categorised to as non-ice-covered lakes (86 lakes).Although some of the seasonally ice-covered lakes may not be completely ice-covered during the cold season and some of the non-ice-covered lakes may have short periods of partial ice cover, the 1 • C lake-mean LSWT offers a good means of evaluating lakes that are typically and non-typically icecovered during the coldest part of the LSWT cycle.To best capture the critical features of both seasonally ice-covered and non-ice-covered lakes, the mean difference in the features between the observed and modelled LSWTs differ with lake category.The tuning approach applied to these two lake categories is summarised in Fig. 2, and described in detail within Sect.2.3.
Using the observed LSWTs (ARC-Lake), the objective of this study is to assess if FLake can be tuned to produce realistic LSWTs for large lakes globally, using relatively few lake properties.It is expected that for each lake, the tuning of lake properties will compensate to a greater or lesser degree for some of the lake-to-lake variability in geographical and physical characteristics.The motivation for this study was to develop a greater understanding of lake dynamics glob-Geosci.Model Dev., 9,2016 www.geosci-model-dev.net/9/2167/2016/ally, offering the potential to help develop parameterisation schemes for lakes in numerical weather prediction models.It is expected that the findings in this study will be of interest to climate modellers, limnologists and current and perspective users of FLake.

Data: ARC-Lake observed LSWTs
LSWT observations from ARC-Lake are used to tune the model.These cover 246 globally distributed large lakes, principally those with surface area > 500 km 2 (Herdendorf, 1982;Lehner and Döll, 2004) but also including 28 globally distributed smaller lakes, the smallest of which is 100 km 2 (Lake Vesijarvi, Finland).The LSWTs are generated from three ATSRs, from August 1991-December 2011(MacCallum and Merchant, 2012).
The ARC-Lake observations have been shown to compare well with in situ LSWT data, demonstrating their suitability for use in this study.Validation of the observations was performed through a match-up data set of in situ temperature data consisting of 52 observation locations covering 18 of the lakes (MacCallum and Merchant, 2012).The timing of ice-on and ice-off events is observed to be consistent with in situ measurements.This is demonstrated through analy-

Model: FLake lake model
FLake is a two-layer parametric representation of the evolving temperature profile of a lake and is based on the net energy budgets (Mironov, 2008).The depth and temperature of the homogeneous "upper mixed layer" and the temperature of the "lake bottom" (representative of the hypolimnion temperature), as illustrated in Fig. 4, are modelled in FLake.The thermal structure of the intermediate stratified layer (thermocline, Fig. 4) is parameterised in FLake through a selfsimilarity representation of the temperature profile, υ(ζ ), using time (T ) and depth (z) as illustrated in Fig. 5.  FLake utilises the minimum set of input data required for 1-D thermal and ice models: meteorological forcing data (shortwave and long-wave radiation, wind speed, air vapour pressure and air temperature), an estimation of turbidity and basic bathymetric data (Lerman et al., 1995).Although models based on the concept of self-similarity are considered to be only fairly accurate (Dutra et al., 2010), we show that mean differences between the model and observed LSWTs are greatly lowered by tuning the model.

Lake-specific model properties
As outlined in the introduction, optimisation of LSWTregulating properties (lake depth (d), snow and ice albedo Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/(α) and light extinction coefficient (κ)) can greatly improve the LSWTs simulated in FLake.Other lake-specific properties adjusted for this study are c_relax_C, fetch, latitude and the starting conditions.c_relax_C: a dimensionless constant used in the relaxation equation for the shape factor with respect to the temperature profile in the thermocline.
The default value of 0.003 was found to be too low to adequately readjust the temperature profile of deep lakes (G.Kirillin, personal communication, 2010), weakening the predicted stratification and affecting the LSWT.For lakes with mean depths < 5 m, the c_relax_C value is set to 10 −2 , and decreases with increasing depth, to a value of 10 −5 for mean depths > 50 m, as recommended by G. Kirillin (personal communication, 2010).
Fetch: wind fetch is calculated as the square root of the product of lake length and breadth for the 205 (of 246) lakes with available dimensions.A line fit on the calculated fetch versus polygon area (in units of km and km 2 respectively) from the GLWD (Lehner and Döll, 2004) of these 205 lakes, showed a strong relationship between fetch and area, Eq. ( 1), R 2 adj = 0.84, p = 0.001.Equation (1), used to determine the fetch of the remaining 41 lakes with no available dimensions, is valid for lakes > 100 km 2 .Although the shape of a lake and its orientation in relation to wind direction are likely to affect wind fetch, this approach is expected to provide reasonable estimates of fetch.fetch = 39.9 km + 0.00781 × area × km −1 (1) Latitude: the latitude of the lake centre reference coordinates (Herdendorf, 1982;Lehner and Döll, 2004).
Starting conditions: provide the initial lake-specific upper mixed layer temperature and depth, bottom temperature, ice thickness and air-ice interface temperature, were shown to shorten the model spin-up time (to an average of < 3 days).A good estimation of the starting conditions was obtained from the FLake model based on the hydrological year 2005/06 (Kirillin et al., 2011).

Fixed model parameters
The icewater_flux, inflow from the catchment and heat flux from sediments remain fixed throughout the investigative and tuning process, across all lakes.For icewater_flux, (heat flow from water to ice) G. Kirillin (personal communication, 2010) suggested values of ∼ 3-5 W m −2 .In this study a value of 5 W m −2 is applied to all lakes.Inflow from the catchment and heat flux from sediments are not considered.

Model forcing data
FLake is forced with ECMWF Interim Re-analysis (ERA) data (Dee et al., 2011;ECMWF, 2009), at the grid points (0.7 • × 0.7 • resolution) closest to the lake centre, shown in the Supplement.The mean daily values of shortwave so- lar downward radiation (SSRD), air temperature and vapour pressure at 2 m, wind speed at 10 m and total cloud cover (TCC), shown in Table 1, are used to force the model.

Tuning method
A range of factors/values for d, α and κ is determined through the model trials (carried out on 21 seasonally icecovered lakes and 14 non-ice-covered lakes, Fig. 6).These lakes broadly represent the range of lake characteristicslake depth, snow and ice albedo and light extinction coefficient -and have available Secchi disk depth data.

Light extinction coefficients for trial lakes
The light extinction coefficient values for the untuned model trial are derived from Secchi disk depth data, κ sd (in units of m −1 ), obtained from the ILEC database (ILEC, 1999).Five methods of relating κ values to Secchi disk depths (Poole and Atkins, 1929;Holmes, 1970;Bukata et al., 1988;Monson, 1992;Armengol et al., 2003) are compared in Fig. 7.These methods cover a range of different water conditions, from coastal turbid waters (Holmes, 1970) and eutrophic water (tested 1 km from a dam in the Sau reservoir, Spain) (Armengol et al., 2003) to a range of North American lakes of different trophic levels (Monson, 1992).
For Secchi disk depths > 10 m, all methods show a reasonably good comparison between Secchi disk depths and κ, Fig. 7. From Secchi disk depths of 10 to 1 m, the range of results between methods becomes increasingly large.Bukata et al. (1998) showed that Eq. ( 2), based on in situ optical measurements from many stations, adequately described Lake Huron, Lake Superior and Lake Ontario, for Secchi disk depths from 2 to 10 m: where S = Secchi disk depth.

ERA data components and description
(3)

Light extinction coefficients for tuning of all lakes
As many lakes do not have available Secchi disk depth data, an alternative approach is used to provide light extinction coefficients in the tuned model trials and for the tuning of Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/all lakes.A range of 10 optical water types, which essentially describe the attenuation process of ocean water and its changes with turbidity (Jerlov, 1976), is applied.These consist of five optical water types for open ocean, type I, IA, IB, II and III; type I being the most transparent and type III being least transparent and five coastal ocean types (1, 3, 5, 7 and 9) (Jerlov, 1976).The spectrum for these 10 ocean water types are divided (in fractions of 0.18, 0.54, 0.28) into bands represented by three wavelengths: 375, 475 and 700 nm respectively.The 10 ocean water types are renamed herein as κ d1 to κ d10 , the values for which are shown in Table 2.

Tuning of lake depth
Lake depth information was obtained from Herdendorf (1982), the ILEC World Lake Database (http:// wldb.ilec.or.jp/),LakeNet (http://www.worldlakes.org/)and Kourzeneva et al. (2012).The mean depth (Z d1 ) is the recommended depth value for FLake.Where only maximum depth is available (nine lakes), the mean depth is calculated using the average maximum-to-mean depth ratio of lakes with known maximum and mean depths.This ratio is 3.5 for seasonally ice-covered lakes and 3.0 for non-ice-covered lakes.In the tuning process, depth factors (outlined in Table 2) are applied to the lake-mean depth.The tuned depth is referred to as the "effective depth".For lakes with no depth information, the effective depth factors are applied to a depth of 5 m.If the resulting effective depth is too shallow, characterised by early LSWT cooling and/or a high summertime LSWTs: July-August-September (JAS) LSWTs, tuning is repeated using a deeper input depth.
For modelling very deep lakes, a "false depth" is recommended (G.Kirillin, personal communication, 2010).For the six lakes with mean depths in excess of 240 m (ranging from 240 m for Lake Kivu to 680 m for Lake Baikal), false depths of 100-200 m are used in the tuning study.

Tuning of snow and ice albedo
FLake uses two categories of albedo for snow (dry snow and melting snow) and for ice (white ice and blue ice).As the snow cover module is not operational in this version of FLake, the snow and ice albedo are set to the same default value in the albedo module; 0.60 for dry snow and white ice and 0.10 for melting snow and blue ice.An empirical formulation is applied in FLake, where the albedo of a frozen lake surface (α lake ) depends on the ice surface temperature (Eq.4) (Rooney and Jones, 2010), accounting somewhat for seasonal changes in albedo (Mironov and Ritter, 2004).By application of this equation, the blue ice (by means of its albedo) has greater influence on the rate of ice melting, in the warming season.
where α w = white ice albedo (0.60), α b = blue ice albedo (0.10), Cα = Ice albedo empirical coefficient (95.6), T 0 = freezing temperature (K) and T p = the surface temperature (K) from the previous time step.The extinction coefficient values used for modelling solar radiation penetration through white ice and blue ice are 17.1 and 8.4 (in units of m −1 ) respectively, and correspond to the top 0.1 m of the ice layer for clear-sky conditions (Launiainen and Cheng, 1998).
The default snow and ice albedo values are referred to as α1 in this study.During the preliminary trials, a higher albedo (than α1) was shown to delay ice-off, substantially improving the timing of early ice-off, compared to observed LSWTs (demonstrated in Fig. 1a).A higher snow and ice albedo causes more of the incoming radiation to be reflected, resulting in a later ice-off.On this basis, we apply three additional higher albedo values (α2 : α4), shown in Table 2, for tuning seasonally ice-covered lakes.Albedo when discussed throughout this study refers to the albedo of snow and ice.The albedo of water (in liquid phase) is maintained at the default value of 0.07 throughout this study.

Wind speed scaling
Scaling of wind speeds is considered during the trials, as most long-term records of wind speed are measured over land (U land ) and are considered to underestimate the wind speed over water (U water ).For adjusting surface wind speeds over land to wind speeds over sea surfaces, Hsu (1988) recommended the scaling shown in Eq. ( 5).For bodies of water with fetch > 16 km, a scaling of 1.2 applied to over-land wind speeds (measured at a height of 10 m) provide reasonable estimates of wind speeds over sea surfaces (Resio et al., 2003).To find a suitable wind speed scaling, the trial work is carried out using the unscaled wind speed (u1), wind speed factored by 1.2 (u2) and wind speed suggested by Hsu (1988)  speed scalings are determined.
where U water = wind speed over water and U land = surface wind speed over land

Summary of the tuning of the LSWT-regulating properties
Table 2 contains a summary of the factors/values for d, α and κ applied in this tuning study.The model is tuned using the optimal combination of LSWT-regulating properties; 80 possible combinations for seasonally ice-covered lakes and 60 possible combinations for non-ice-covered lakes.

Tuning metrics
The tuning metrics are the mean differences (between the modelled and the observed LSWTs), used to quantify the effect that the LSWT-regulating properties have on the modelled LSWTs.The metrics (normalised and equally weighted) determine the optimal LSWT model selected for each lake.

Tuning metrics for seasonally ice-covered lakes
The metrics and the effect of the LSWT-regulating properties on them, for seasonally ice-covered lakes is summarised in Table 3.The effect of light extinction coefficient on the JAS LSWTs is demonstrated in Fig. 8, showing that the tuned light extinction coefficient (κ d ) value, κ d6 in place of a lower (more transparent) κ d value (κ d2 ), described in Table 2, substantially improves the JAS LSWTs.In this figure, the greater effect of light extinction coefficient on the maximum LSWT

Tuning metrics for non-ice-covered lakes
For these lakes, the daily MAD and the difference between the observations and model for the months, where the minimum and maximum observed LSWTs occur (mth min and mth max ), are applied as metrics.Although there no definitive stages in the LSWT cycle for non-ice-covered lakes, the mth min and mth max exert some control over temporally rec-   showing that when modelled at a greater depth, the lake cools later and the maximum LSWT is lower.
onciling the modelled monthly extremes with the observed monthly extremes.

Additional metrics for seasonally ice-covered lakes and non-ice-covered lakes
For each lake, the fraction of the observed mean LSWT variance over the number of years with observations, accounted for in the tuned model is used to help independently evaluate the tuned LSWTs.For non-ice-covered lakes, the variance is determined using var min (and var max ): the variance in the mean LSWT for the month in which the minimum (and maximum) LSWT is observed.For seasonally ice-covered lakes, the variance is determined using var jas : the variance in the observed mean JAS LSWTs.The fraction of these observed LSWT variances accounted for in the tuned model are quantified; inter min , inter max and inter jas (using R 2 adj ).The calculations to quantify var jas and inter jas are shown in Eqs. ( 6) to (8).var jas : observed JAS LSWT variance (in units of K 2 ) over the length of the tuning period (Walker and Shostak, 2010): where x obs−jas = observed mean JAS LSWTs for each individual year x = mean across all years, N = number of years with JAS LSWTs inter jas : the fraction of the observed JAS LSWT inter-annual variance (var jas ) accounted for in the tuned model (Lane et al., 2016): where P = total number of regressors (Walker and Shostak, 2010) where obs − jas = observed JAS LSWT and mod − jas = modelled JAS LSWT.
The same equations, Eqs. ( 6) to (8), are applied to determine Inter max , var max , inter min and var min , substituting "JAS" with "max" and "min", where obs − min (and modmin) = mean-observed LSWTs (and modelled LSWTs) in the month where the minimum LSWT occurs, and where www.geosci-model-dev.net/9/2167/2016/Geosci.Model Dev., 9, 2167-2189, 2016 Table 4.The effect of wind speed scaling on untuned modelled LSWTs, presented as the mean difference, between the modelled and observed values, across lakes with the spread of differences defined as 2σ , where wind speeds u1 is unscaled, u2 is factored by 1.2 and u3 (U water = 1.62 m s −1 + 1.17U land ).Results are presented for seasonally ice-covered and non-ice-covered trial lakes.Results highlight that u3 is most applicable to seasonally ice-covered lakes but there is no one wind speed most suited for all lakes (while the mean difference is improved with u3, the spread of the mean differences across lakes for mth min and mth max show little change).Units are given in square brackets.
Trial results for untuned model obs-max (and mod-max) = mean-observed LSWTs (and modelled LSWTs) in the month where the maximum LSWT occurs.
3 Trial results for wind speed scaling For both seasonally ice-covered lakes and non-ice-covered lakes, wind speeds, u1, u2 and u3 were modelled with untuned LSWT properties: mean lake depth (Z d1 ), default snow and ice albedo (α1) and light extinction coefficient derived from Secchi disk depth data (κ sd ).The trials show that wind speed has a consistent effect on the modelled LSWT of seasonally ice-covered lakes.The higher wind speed scaling (u3) causes earlier cooling and later warming (reducing the 1 • C cooling day and 1 • C warming day mean differences), lengthening the ice cover period and lowering the JAS LSWT, as demonstrated for Lake Simcoe, Canada, in Fig. 10.It is expected that the tuning of d, α and κ, with an applied wind speed of u3, will produce modelled LSWTs substantially closer to the observed LSWTs than those shown in Fig. 10, where tuning of d, α and κ is not applied.The more rapid mixing and heat exchange between the surface and atmosphere, as a result of the higher wind speed, causes an earlier modelled 1 • C cooling day.As wind promotes ice growth in the model, higher wind speeds also contribute to the later-modelled 1 • C warming day.Wind speed scaling, u3 in place of u1, for the trial seasonally ice-covered lakes, reduces the mean difference in the length of the average cold phase (when compared to the observed cold phase) by ∼ 50 % (from 39 to 21 days) and reduces the JAS LSWT mean difference by ∼ 50 %, from 3.71 to 1.87 • C, Table 4.
For non-ice-covered trial lakes, five of the seven lakes at latitudes > 35 • N/S (north or south) show best results with u3, as demonstrated for Lake Biwa, located at 35.6 • N, Fig. 11a.Five of the seven lakes located < 35 • N/S show the best results with u1, as demonstrated for Lake Turkana, located at 3.5 • N, Fig. 11b.
For the remainder of the tuning trials (tuned using the range of d, α and κ factors/values outlined in Table 2), wind speed scaling u3 was applied to all seasonally ice-covered and to non-ice-covered lakes > 35 • N/S and no scaling (u1) to non-ice-covered lakes < 35 • N/S.As the daily MAD target of < 1.0 • C is met for both lake categories (shown in the second results column in Table 5), the tuning approach described here is applied to all lakes.

Summary of results
The daily MAD and spread of differences (2σ ) between the modelled and observed LSWTs, across the seasonally ice-covered lakes and non-ice-covered lakes, is reduced from 3.07  5.The tuned values for the LSWTregulating properties for all lakes and the tuning metrics are shown in the Supplement.
These results demonstrate that the tuning process with the applied wind speed scalings can provide significant improvements on the untuned model: run using the lake mean depth, light extinction coefficients derived from Secchi disk depth (as shown in Sect.2.3.1) and the model default albedo (seasonally ice-covered lakes only).
The applied tuning method yielded a daily MAD of 0.74 • C ± 0.48 • C, across 135 of the 160 seasonally icecovered lakes, Table 6.The remaining 25 seasonally icecovered lakes yielded comparatively poor results; the 1 • C cooling day was 14 days too early and/or the JAS LSWT was ≥ 2 • C higher than observed LSWTs.Re-tuned using greater effective depth factors and higher κ d values, as outlined in the next sub-section (Sect.4.1.1),yielded a daily MAD of 1.11 • C ± 0.56 • C, across the 25 lakes.Across the 160 lakes, a daily MAD of below 1 • C was achieved (0.80 • C ± 0.56 • C, Table 5).A daily MAD of below 1 • C is again achieved (0.96 • C ± 0.66 • C) when 84 of the 86 non-ice-covered lakes are considered (Table 5).The remaining two lakes yielded highly unsatisfactory results.

Seasonally ice-covered lakes
The 25 lakes that were re-tuned are shallow (average mean depth < 5 m) and small (18 of the 25 lakes are < 800 km 2 ), relative to the depth and area of the larger seasonally icecovered lakes.Twenty (20) of the 25 lakes are located in eastern Europe or Asia, at relatively low altitudes; 22 of the 25 lakes are < 752 m a.s.l.On initial tuning, these 25 lakes were tuned to the highest depth factor, Z d4 (1.5 times the mean depth) and/or the highest light extinction coefficient, Figure 11.Effect of wind speed scaling on lake surface water temperatures (LSWTs) for a temperate non-ice-covered lake (a) Lake Biwa, Japan (36 • N 136 • E), and for a tropical non-ice-covered lake (b) Lake Turkana, Africa (4 • N 36 • E), showing that the modelled LSWT for the temperate lake is better represented using u3 (U water = 1.62 m s −1 + 1.17U land ), and the modelled LSWT for the tropical lake is better represented using u1 (unscaled wind speed).mth min (and mth max ) is the difference between the observed and modelled LSWTs for the month where the minimum (and maximum) LSWT is observed.κ d5 (lowest transparency).Although the transparencies for these 25 lakes are largely unknown, shallow lakes generally have poorer light transparencies than deeper lakes due to upwelling of bottom sediment.The shallow depth of the modelled lake (lower heat capacity) and the poor transparency of water (more heat retained in surface) were evident in the metric results: early 1 • C cooling day and/or high JAS LSWT values.A modified tuning set-up, to allow for a greater modelled depth to increase the heat capacity -postponing the 1 • C cooling day -and lower transparency values (higher κ d ), causing less heat to be retained in the surface and lowering the JAS LSWT, is applied to these 25 lakes.
The tuning set-up modified to include three greater depth factors of 2.5, 2 and 4 times the mean depth (Z d6 , Z d7 and Z d8 ) and two higher light extinction coefficient values, κ d6 and κ d7 (Table 2), substantially improves the 1 • C cooling day and the JAS LSWT for these 25 lakes.A summary of the results are shown in Table 6 column 2. The tuning metric results for all 160 lakes (using the modified tuning set-up for the 25 shallow lakes) are illustrated in Fig. 12.
www.geosci-model-dev.net/9/2167/2016/Geosci.Model Dev., 9, 2167-2189, 2016 Table 5. Summary of the untuned and tuned metrics for the trial lakes and the tuned metrics for all lakes (metrics are explained in Table 4).The results, presented for seasonally ice-covered and non-ice-covered lakes in each instance, show the mean between the modelled and observed values, across lakes with the spread of differences defined as 2σ .For tuned lakes, wind speed scaling u3 was applied to all seasonally ice-covered and to non-ice-covered lakes > 35

Non-ice-covered lakes
The tuning metric results for each of the 84 lakes are illustrated in Fig. 13 and a summary of these results are shown in Table 5.The poor tuning results, observed for 2 of the 86 lakes (Lake Viedma and the Dead Sea) are most likely due to differences between the altitude of the ERA T 2 air temperature (geopotential height) and the lake altitude.Lake Viedma, an Argentinian freshwater lake of unknown depth, yielded a daily MAD of 3.1 • C. The Dead Sea, a deep and highly saline (340 g L −1 ) lake located in Asia at 404 m b.s.l., yielded a daily MAD of 4.1 • C.
For the Dead Sea, a temperature difference (in the month of maximum temperature) between the observed LSWT (33 • C) and ERA T 2 air temperature (25 • C), results in a negative modelled mean difference of 6.3 • C in LSWT for this month.Given the standard air temperature lapse rate (6.5 • C km −1 ), altitude can explain the substantially lower air temperatures.The altitude of Dead Sea (−404 m a.s.l.) is lower by ∼ 850 m a.s.l.than the altitude of the meteorological data at the lake centre co-ordinates, 445 m a.s.l.(determined by interpolating surrounding cells using the orography data accompanying the ECMWF meteorological data).
For Lake Viedma, while the observed LSWTs range from 5 to 10 • C, the minimum ERA T 2 air temperature remains well below 0 • C for many months of the year, regularly reaching −8 • C, resulting in a negative modelled mean difference of 4.8 • C for the month of minimum LSWT.This difference can be, at least, partially explained by the difference in altitude (> 500 m a.s.l.) between the altitude of Lake Viedma (297 m a.s.l.) and the altitude of meteorological data (825 m a.s.l.) at the lake centre co-ordinates.

Tuning of saline and high altitude lakes
The tuned metrics shown in Table 7 (seasonally ice-covered lakes) and in Table 8 (non-ice-covered lakes) indicate that FLake is successful for tuning both saline and high altitude lakes, as well as freshwater and low altitude lakes.
Although the density of freshwater in FLake is determined at sea level (normal atmospheric pressure) (Mironov, 2008) and the altitude of lakes are not directly considered in FLake, lake altitude (ranging from −12 to 5000 m a.s.l., over the 246 lakes) is considered indirectly through the altitude of Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/   the meteorological forcing data (ERA) at the lake centre coordinates.
The majority of the high altitude lakes are also saline: 7 of the 10 non-ice-covered lakes and 12 of the 14 seasonally ice-covered lakes.The good agreement between the observed and modelled LSWTs for two high altitude lakes (> 1500 m a.s.l.) is shown in Fig. 14.

Independent evaluation
Two methods are used to independently evaluate the tuned model.
The fraction of observed LSWT variance that is detected in the tuned model is quantified (using R 2 adj ); inter min and inter max (non-ice-covered lakes) quantifies the observed variance in the month in which the minimum LSWT (var min ) and maximum LSWT (var max ) occurs and inter jas (seasonally ice-covered lakes) quantifies the observed variance in the mean JAS LSWT (var jas ).
The metrics for 2011 (not used in the tuning process) are compared with metrics from 2 tuned years.

Variance detected in the tuned model
The results show that the modelled LSWTs capture less of the true (observed) inter-annual variance in lakes where the observed LSWT variance and the annual LSWT range is low.This indicates that lower latitude lakes and high altitude lakes are less well simulated in the model, than lakes with greater observed LSWT variance and annual range.This would also indicate that lakes in the Southern Hemisphere at 35-55 • S are less well simulated than lakes in the Northern Hemisphere at the same latitude, as the annual LSWT range is considerably lower at 35-55 • S than at 35-55 • N (Layden et  al., 2015).
For non-ice-covered temperate lakes, the inter max and inter min fraction is substantially greater (0.49 and 0.37) than in tropical lakes (0.07 and 0.13), Table 9.This can be ex-Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/ports the findings that tropical and high altitude lakes are less well simulated in the model.

Comparison of tuned and untuned model LSWTs
The final year (2011) of available observational ARC-Lake LSWT data is used to independently evaluate the tuning process.The metrics from the tuned model forced for the year 2011 are compared with the metrics from 2 tuned years (1996; the first full year of data from ATSR2 and 2010; the last year of tuned data from Advanced ATSR (AATSR), as shown in Tables 10 and 11.The mean metric results and the spread of differences across the 135 seasonally ice-covered lakes are highly comparable across all 3 years, with marginally better daily MAD metrics observed for the untuned year.For the 25 shallow lakes tuned with the modified tuning set-up, the MAD result for the untuned year is comparable with 2010.For the other three metrics, the untuned year has a lower spread of differences across lakes than those for 2010.The 1 • C warming day difference for the untuned year is greater than the difference for 2010 but less for 1996.The 1 • C cooling and warming day mean differences for 1996 and 2010 are less comparable for the 25 lakes than for the 135 lakes.This may be because the modelled effect of depth on the metrics is more predictable for deeper lakes, as illustrated in Fig. 15, than for shallow lakes.
Although inter-annual variance may somewhat obscure year-on-year comparisons, the results of the modelled LSWTs for the untuned year (2011) compare well to the modelled results from the tuned years (1996 and 2010) showing that the model remains stable when run with ERA forc- ing data outside the tuning period.For non-ice-covered lakes, although the daily MAD and dispersion of errors is slightly higher for the untuned year (2011, Table 11), overall the metrics are very comparable to the metrics from 1996 and 2010.
5 Findings and discussion 5.1 The effect of the 1 • C warming day on JAS LSWT Through the trial work, the modelled effect of a change in snow and ice albedo on the timing of 1 • C warming day (indicative of ice-off), the JAS LSWT and the timing of the 1 • C cooling day (indicative of ice-on) is demonstrated, for deep high latitude or very deep seasonally ice-covered lakes.
Using the default snow and ice albedo (α1, Table 2), the modelled 1 • C warming day of the 21 trial lakes occur, on average, 20 days too early.Across these lakes, a higher snow and ice albedo (α2, Table 2) results in a delay to LSWT warming to 1 • C by 27 days ± 12.6 days and a decrease in the JAS LSWT difference (between modelled and observed LSWTs) by ∼ 50 %, to 0.98 • C ± 2.51 • C. Much of the modelled variance in the JAS LSWT decrease, across the 21 lakes, is attributed to lake depth and latitude; accounting for 0.50 (R 2 adj , p = 0.001) of the variance (stepwise regression).Separately, depth accounts for 0.35 (p = 0.003) and latitude for 0.26 (p = 0.01) of the variance.
The LSWTs for two deep high latitude lakes (Great Bear and Great Slave) modelled with α2 (high) and α1 (low; default) snow and ice albedos, shown in Fig. 16, clearly shows the modelled effect of snow and ice albedo on the warming day and on the JAS LSWT.While the modelled change in the snow and ice albedo is the cause of the delay in the 1 • C warming day, we find that the decrease in the JAS LSWT for Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/Table 10.Results of independent evaluation of the tuning process for seasonally ice-covered lakes.The spread of differences across lakes is defined as 2σ .These results illustrate that the metrics (explained in Table 4) from the untuned year (2011) compare well with metrics from 1996 (the first full year of data from along-track scanning radiometers 2 (ATSR2)) and 2010 (the last year of tuned data from Advanced ATSR.For the untuned year ( 2011)) for each lake, the model is forced with the effective lake depth (Zd), snow and ice albedo (α) and light extinction coefficient (κ d ) values determined during the tuning process, shown in the supplement.Units are given in square brackets.Table 11.Results of the independent evaluation of the tuning process for non-ice-covered lakes.The spread of differences across lakes is defined as 2σ .Metrics (explained in Table 4) for the untuned year (2011) are compared with those from the first full year of data from along-track scanning radiometers 2 (ATSR2) (1996) and the last year of tuned data from Advanced ATSR (2010).For the untuned year (2011), for each lake, the model is forced with the effective lake depth (Z d ), snow and ice albedo (α) and light extinction coefficient (κ d ) values determined during the tuning process, shown in the supplement.Units are given in square brackets.A study on Lake Superior (average depth of 147 m) shows a JAS LSWT warming trend (of 2.5 • C from 1979 to 2006) substantially in excess of the air temperature warming trend (Austin and Colman, 2007).Austin and Colman attributed this warming trend to a longer warming period, caused by an earlier ice-off date, of ∼ 0.5 day yr −1 .Foster and Heidinger (2014) suggested warming trends in North America may be due to changes in cloud albedo; with an observed loss of 4.2 % in total cloudiness between 1982 and 2012.
As shown, the modelled decrease in the JAS LSWT (as a result of the higher albedo; α2) is more pronounced for deep lakes.The modelled 1 • C cooling day is also shown to occur earlier in these lakes, with deeper lakes showing a greater earlier cooling.The 1 • C cooling day occurs 3.4 days earlier for Great Slave, average depth of 41 m (decrease of 4.26 • C in the modelled JAS LSWT).For Great Bear, average depth of 72 m, which shows a modelled JAS LSWT decrease of 3.40  138 m, shows a modelled JAS LSWT decrease of 2.60 • C and an earlier 1 • C cooling day, by 12.8 days.
The findings are sensible.A delay in the 1 • C warming day, shortening the lake warming period, may not prevent a shallow lake reaching its full heating capacity but may prevent a deep lake from reaching its maximum heat storage capacity.At higher latitudes, the LSWT warming period for northern hemispheric lakes becomes increasingly short (Layden et al., 2015).As a result, deep lakes increasingly fall short of reaching their maximum heat storage, causing a larger JAS LSWT decrease.Any changes to the 1 • C warming day of deep and high latitude (or high altitude) lakes will therefore affect JAS LSWT.Deep lakes also cool more slowly than shallow lakes, resulting in a later cooling day.
These findings highlight the sensitivity of the whole LSWT cycle of deep high latitude lakes, to changes in snow and ice albedo and in the timing of the 1 • C warming day, as illustrated in Fig. 15.This figure also illustrates how an earlier 1 • C cooling day caused by a lower JAS LSWT may be counteracted or masked in deep lakes, where heat is retained during the cooling period.
The effect that depth has on the JAS LSWT is apparent when comparing lakes at the same altitude and latitude but with different depths.For example, Lake Nipigion and Lake Manitoba, both located in Canada (50 and 51 • N) and at similar altitudes (283 and 247 m a.s.l) have considerably differ-ent depths, 55 and 12 m respectively.Significant differences are observed in JAS LSWT for these lakes, the deeper lake having an average JAS LSWT 4.4 • C lower than that of the shallower lake (15.4 • C compared to 19.8 • C).
As the snow cover module with FLake is not operational in this version of the model; the insulating effect that snow has on the underlying ice is not modelled.This may contribute to the earlier 1 • C warming day and the higher JAS LSWTs, when modelled with the default albedo.It is also possible that the icewater_flux value of 5 W m −2 may be an overestimation of the water-to-ice heat flux in the ice growth phase of deep and shallow lakes.This greater heat flux, leading to underestimated ice thickness, could have contributed to the large 1 • C warming day mean difference shown in Table 5 (column 1).In a study by Malm et al. (1997), the water-toice heat flux during the ice growth phase was shown to be < 1 W m −2 in both deep (15-20 m) and shallow lakes.Underestimated ice thickness, causing an early ice melt, may possibly have led to over-tuning of albedo in the tuned model.

Lake-bottom temperatures modelled in FLake
The month of minimum LSWTs in the annual cycle (monthly minimum) have the potential to be used as a proxy for determining the temperature of the bottom layer (hypolimnion) of non-ice-covered lakes.The monthly minimum climatological ARC-Lake LSWT explains 0.97 of the inter-lake variance (using R 2 adj ) in the bottom temperatures, obtained from the FLake model based on the hydrological year 2005/2006 (Kirillin et al., 2011) and have a ∼ 1 : 1 relationship, as shown in Fig. 18.Although FLake is a two-layer model, the depth of the hypolimnion layer is not calculated; the bottom modelled temperature is representative of the hypolimnion temperature, which remains constant with depth.
Empirically, it has previously been shown that from the Equator to approximately 40 • (N/S), the steep decline in the minimum LSWT is reflected in the hypolimnion temperature (Lewis Jr., 1996).This relationship is applicable to deep stratified non-ice-covered lakes.For these lakes, the surface water, when at its coolest in the annual cycle (minimum LSWT) and therefore its densest, sinks to the lake-bottom.During the summer stratification period, the water in the upper mixed layer is warmer and less dense and therefore remains in the upper layer (with exception to high wind or storm conditions, which can induce intense vertical mixing).The strengthened density gradient in the summer thermocline (as demonstrated for Lake Malawi in Fig. 4) also protects the hypolimnion from heat flux through the lake surface.As a result, the lake hypolimnion temperature of deep non-ice-covered lakes can reflect the minimum LSWT.The good agreement between the monthly minimum LSWT (using the ARC-Lake monthly minimum climatology LSWTs) and the bottom temperature, for all deep (> 25 m) non-icecovered lakes (14 lakes) supports this empirical observation (Fig. 18).
Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/In FLake, the bottom temperature is not independent of surface temperature; the change in the surface heat flux over time is used in calculating the upper mixed layer temperature, and the difference in heat flux between the upper mixed layer and lake bottom are considered in the lake-bottom temperature calculation (Kourzeneva and Braslavsky, 2005).Although the minimum surface temperature is therefore related to the bottom temperature in FLake, the good agreement between minimum ARC-Lake LSWTs and the bottom temperatures indicate that the monthly minimum LSWTs are a potential proxy for determining the lake-bottom temperature.This also supports Lewis's empirical relationship between lake surface temperature and lake-bottom temperature.
Although changes in other factors affect hypolimnion temperature, such as influx of cooler water and geothermal heating, the monthly minimum LSWTs from satellites can offer a good indication of hypolimnion temperature; useful in cases where this otherwise cannot be or is not observed directly.

Wind speed scaling for low latitude lakes
There is no optimal scaling for all non-ice-covered lakes, as discussed in Sect.3.This is possibly attributable to the highly variable range of latitudes, LSWTs and mixing regimes of non-ice-covered lakes.
For the deep (> 25 m) non-ice-covered lakes (14 lakes), the density difference between the lake surface (in the month of maximum LSWT) and the hypolimnion during the summer stratification period when the density (and temperature) gradient of the thermocline is strongest, as illustrated in Fig. 4, was calculated (Haynes, 2013).For lakes at latitudes below 35 • N/S, the average density difference between these two layers is substantially lower (0.352 × 10 −3 kg m −3 ) than for lakes at latitudes above 35 • N/S (1.183 × 10 −3 kg m −3 ).This is due to the smaller annual temperature range of the lower latitude lakes.
It is possible that the large density difference between the lake surface at maximum LSWT and the hypolimnion in high latitude lakes during the stratification period may produce a buffer against wind-induced mixing and therefore lessen the heat flux through the thermocline.As winds can drive lake mixing in deep lakes, it strongly influences the mixed layer depth and the LSWT.The larger the temperature (and density) gradient between the lake surface and the hypolimnion during stratification, the more wind energy is required to produce the same amount of mixing than for lakes with a smaller temperature (and density) gradient between the two layers.Although the density differences between the two layers are considered in FLake, the model is forced with over-land wind speed measurements.It is possible that when forced with an underestimated wind speed, the effect of wind on the LSWT will be further reduced.As a result, higher latitude lakes may show more representative LSWTs using a higher wind speed scaling.

Improving modelled LSWTs in FLake
The optimal LSWT-regulating properties of the 244 lakes provide a guide to improving the LSWT modelling in FLake for other lakes, without having to tune the model for each lake separately.

Depth
The tuning results show that deep lakes are generally tuned to a shallower effective depth and shallower lakes to a deeper effective depth.Figure 19 shows the relationship between the lake-mean depth and the effective (tuned) depth of all 244 successfully tuned lakes, colour coded by the effective depth factor optimised in the tuning process.The figure legend shows that the effective depth factor decreases with increasing average lake depth (also graphed in the figure insert), providing a means to estimate an appropriate effective depth for any lake with a mean depth from 4 to 124 m.
The tuned lake depths are sensible.For shallow lakes, tuning to a deeper effective depth may compensate for not having considered the "heat flux from sediments" scheme in the model.Retention of heat in the sediments of a lake has the same effect on modelled heat storage capacity, as deepening the lake.
Many deep lakes have three distinct layers, the upper mixed layer, the underlying thermocline and the bottom layer (hypolimnion), illustrated in Fig. 4. As FLake is essentially a two-layer model, it is possible that for deep lakes the mean depth (mean of entire lake depth) is tuned to a shallower effective depth as it is more representative of the mean depth of the two upper lake layers.Other factors affecting the rate at which heat is exchanged between the atmosphere and the www.geosci-model-dev.net/9/2167/2016/Geosci.Model Dev., 9, 2167-2189, 2016 surface water, such as topography, altitude, bathymetry and surface area are not considered in FLake.It is possible that lake depth tuning may also compensate for the effect that these factors have on the rate of the surface heat exchange.

Light extinction coefficient
Across all lakes, 57 % were tuned to light extinction coefficient values of κ d4 or κ d5 .These lakes are globally distributed and have a wide range of mean depths (1-138 m) with an average mean depth of 16 m.In view of this finding and considering that light extinction coefficient values are scarce for the majority of lakes, we assess if κ d4 and κ d5 can be used to provide a good estimation of the light extinction coefficient for modelling LSWTs in FLake.
The untuned model is forced using two sets of light extinction coefficient values and the daily MAD results are compared.In the first model run, the average κ sd value (derived from Secchi disk depth data) of the trial lakes of each lake type is applied to all lakes of corresponding type: for the 21 seasonally ice-covered trial lakes, κ sd = 0.82 m −1 ; for the 14 non-ice-covered trial lakes, κ sd = 1.46 m −1 .In the second run, the model is forced with κ d4 or κ d5 values: κ d4 is applied to all lakes > 16 m in depth (the average depth of lakes tuned with κ d4 or κ d5 ) and κ d5 to all lakes < 16 m in depth.It makes practical sense to apply the less transparent of these two κ d values (κ d5 ) to shallower lakes, as shallow lakes are generally more affected by lake-bottom sediments than deeper lakes.
For both model runs the default albedo and the mean depth are applied, while all other model parameters are kept the same.A comparison of the two model runs shows that when LSWTs are modelled with κ d4 and κ d5 values, the daily MAD is reduced from 3.38 • C ± 2.74 • C to 2.28 • C ± 2.30 • C (33 % decrease) across all lakes.This in-dicates that in the absence of available light extinction coefficient values, application of κ d4 and κ d5 values may improve the modelling of LSWTs of large lakes in FLake.

Snow and ice albedo
For seasonally ice-covered lakes, only 19 % of the lakes were tuned to the default snow and ice albedo, α1, (snow and white ice = 0.60 and melting snow and blue ice = 0.10); 64 % of lakes were tuned to two higher albedos α2 or α3, (snow and white ice = 0.80 and melting snow and blue ice = 0.60 for α2 or 0.40 for α3), indicating that the default snow and ice albedo may be too low for the majority of lakes.In the absence of lake-specific snow and ice albedo information, the albedo value α3 (snow and white ice = 0.80, melting snow and blue ice = 0.40) may provide a good estimate.The α3 values are highly comparable to albedo values measured on a lake in Minnesota using radiation sensors, where the mean albedo of new snow was shown to be 0.83 and the mean ice albedo (after snow melt) was 0.38 (Henneman and Stefan, 1999).

Summary and conclusions
The 1-D freshwater lake model, FLake, was successfully tuned for 244 globally distributed large lakes (including saline and high altitude lakes) using observed LSWTs (ARC-Lake), for the period 1991 to 2010.This process substantially improves the measured mean differences in various features of the lake annual cycle, using only three lake properties (depth, snow and ice albedo and light extinction coefficient), as summarised in Table 5.In the process of tuning the model, we demonstrate several aspects of LSWT behaviour, in a way that cannot be done using the LSWT observations alone.We demonstrate the dependency of the whole modelled LSWT cycle of deep high latitude or high altitude lakes, on changes in the snow and ice albedo.The monthly minimum LSWTs from satellites are demonstrated to offer a good indication of the modelled lake-bottom temperature, with a 1 : 1 relationship shown (Fig. 18).This is highly useful where the lakebottom temperature cannot be or is not observed directly.
The amount of observed inter-annual LSWT variance (in the month in which the minimum LSWT and maximum LSWT occurs for non-ice-covered lakes and in the JAS LSWTs for seasonally ice-covered lakes) detected in the tuned model was quantified.It can be concluded that lakes where the observed LSWT variance is low (lower latitude and high altitude) and for non-ice covered where the annual range is low are less well simulated in the model, than lakes with greater observed LSWT variance and annual range.
We found that no wind speed scaling, u1, is most appropriate for lakes at lower latitudes, < 35 • N/S, and that the largest scaling (u3; U water = 1.62 m s −1 + 1.17U land ) is most appropriate for lakes at higher latitudes > 35 • N/S.A greater Geosci.Model Dev., 9, 2167-2189, 2016 www.geosci-model-dev.net/9/2167/2016/resistance to wind-induced mixing and heat flux through the thermocline, as a result of a greater density gradient between the lake surface and the hypolimnion of high latitude lakes, may explain the suitability of the largest scaling to higher latitude lakes.The optimal LSWT-regulating properties (effective depth, snow and ice albedo and light extinction) for the 244 lakes are shown to be sensible and may provide a guide to improving the LSWT modelling in FLake for other lakes, without having to apply a tuning process to the model.
The relationship between the lake-mean depth and the effective (tuned) depth of all 244 successfully tuned lakes, show that deep lakes are generally tuned to a lower depth and shallower lakes to a greater depth.Figure 19 provides a means to estimate an appropriate effective depth for any lake with a mean depth from 4 to 124 m.An albedo value α3 (snow and white ice = 0.80, melting snow and blue ice = 0.40) is recommended in place of the default value (α1).Where κ values are unknown, applying κ d4 for lakes > 16 m in depth and κ d5 for lakes < 16 m in depth improves the modelled LSWT.
This paper predominantly focused on the tuning of FLake and interpretation of the LSWT annual cycle using the tuned model.The tuned model is forced with ERA data over the available time span of LSWT observations (16-20 years) but has the potential to be forced with ERA data covering a longer time span (ERA data are available for a period of > 33 years; 1979-2012).This offers the potential to provide a better representation of LSWT changes over a longer period of time, as satellite observations for the relatively short period may reflect some inter-annual variance.As demonstrated, the use of remote sensing and modelled LSWTs together extend the reliable quantitative details of lake behaviour beyond the information from either remote sensing or models alone.The ARC-Lake data set has since been extended to include ∼ 1000 smaller lakes (surface area > 100 km 2 ) worldwide, offering the potential to further quantify aspects of lake behaviour worldwide.
The findings in this study are expected to be of interest to limnologists concerned with the relationship between certain features of the LSWT cycle and lake characteristics.Limnologists may also benefit from other aspects of this study, for example, the effect of wind speed scaling on LSWTs and how the observed minimum monthly LSWTs may be used to estimate lake-bottom temperatures.The optimal LSWTregulating properties of the 244 lakes may provide a guide to current and prospective users of FLake for improving the LSWT modelling in FLake for other lakes, without having to tune the model for each lake separately.This is of particular use for lakes where lake characteristic information is not available.The described approach to this study can provide practical guidance to scientists wishing to tune FLake to produce reliable LSWTs for new lakes.

Figure 1 .
Figure1.Preliminary modelled runs for Lake Athabasca, Canada (59 • N 110 • W), showing that adjustments to lake depth (d), snow and ice albedo (α) and light extinction coefficient (κ) can greatly improve the modelled lake surface water temperatures (LSWTs) compared to the default/recommended d, α and κ values; (a) shows that a higher α causes a later ice-off date, comparing well with the observed (ARC-Lake) ice-off date, (b) shows that a lower d causes an earlier ice-on date and a lower κ value (greater transparency) reduces the maximum LSWT and (c) shows that the combined effect of the adjusted d, α and κ produce LSWTs that are highly comparable to the observed ARC-Lake LSWTs.

Figure 2 .
Figure2.Study approach overview (trials, tuning, evaluation and results) for (a) seasonally ice-covered lakes and (b) non-ice-covered lakes.For the trials, wind speed scaling, u1, u2 (recommended for lakes with fetch > 16 km) and u3 (recommended for open ocean water) is assessed on the untuned model, tuning is then trialled with a range of factors for d and values for α and κ using the selected wind speed scaling.The tuning approach produces modelled LSWTs for all possible combination of d, α and κ, 80 modelled runs for seasonally ice-covered lakes and 60 for non-ice-covered lakes.For the evaluation, the tuning metrics (normalised and equally weighted) are the basis for selection of the optimal (tuned) LSWT model for each lake.

Figure 3 .
Figure 3. Location of 246 observed lakes colour coded by surface area (obtained using polygon area in Global Lakes and Wetlands Database (GLWD) showing zoomed inset of North America and northern Europe.

Figure 4 .
Figure 4. Winter and summer depth and temperature profile for Lake Malawi (mean depth of 273 m), Africa (12• S 35 • E), illustrated using data from the ILEC world lake database (http://wldb.ilec.or.jp/), showing the three distinct layers (in winter and summer) of a deep stratified non-ice-covered lake.FLake, a two-layer model, predicts the depth and temperature of the "upper mixed layer" and the temperature of the "bottom layer" (shown on the left), and "thermocline" depth and temperature profile (shown on the right).

Figure 5 Figure 5 .
Figure 5 Schematic representation of the temperature profile in the upper mixed layer and in the thermocline, reproduced from Killirin (2003).The self-similarity representation of the temperature profile in FLake is determined using dimensionless co-ordinates, ζ = (zh)/Δh, and υ = [T(z)-T D ]/ΔT.

Figure 6 Figure 6 .
Figure 6 Location of lakes, with red square showing the trial lakes a) 160 seasonally ice-covered lakes, including 21 trial lakes and b) 86 non-ice covered lakes including 14 trial lakes

Figure 7 .
Figure 7.A comparison of five methods relating light extinction coefficients to Secchi disk depths, showing that all methods compare reasonably well at Secchi disk depths > 10 m.

)
Figure9.Effect of depth on the lake surface water temperature (LSWT) for Lake Ladoga, Russia (61 • N 31 • E) (mean depth 52 m), showing that when modelled at a greater depth, the lake cools later and the maximum LSWT is lower.

Figure 10 .
Figure 10.Effect of wind speed scalings on the modelled lake surface water temperature (LSWT) for Lake Simcoe, Canada, 44 • N 79 • W (depth 25 m), showing that the greatest wind speed scaling, u3 (U water = 1.62 m s −1 + 1.17U land ), in place of the unscaled wind speed, u1, reduces the daily mean absolute difference and July-August-September LSWT mean difference by ∼ 50 %.Modelled with untuned LSWT properties: mean lake depth (Z d1 ), default snow and ice albedo (α1) and light extinction coefficient are derived from Secchi disk depth data (κ sd )

Figure 12 .
Figure 12.Tuning metric mean differences between modelled and observed LSWTs for all 160 lakes with seasonal ice-cover.The results for the 25 lakes tuned with modified tuning approach are marked by diamond symbols (a) July-August-September (JAS) LSWT mean difference, (b) daily mean absolute difference (MAD), (c) 1 • C cooling day mean difference and (d) 1 • C warming day mean difference.

Figure 13 .
Figure 13.Tuning metric results for the 84 non-ice-covered lakes (a) daily mean absolute difference (MAD) between observed and modelled LSWTs, (b) mth max and (c) mth min .mth min (and mth max ) is the difference between the observed and modelled LSWTs for the month where the minimum (and maximum) LSWT is observed.

Figure 15 Figure 15 .
Figure 15 Schematic linking the modelled interactions between the lake surface water temperature (LSWT) regulating parameters: lake depth (d), snow and ice albedo (α) and light extinction coefficient (κ), shown in squares and wind (shown in triangle) with the LSWT metrics: 1 º C cooling day, 1 º C warming days and July August September (JAS) LSWT (shown in circles).
± 0.81 mth max [ • C] −0.23 ± 2.40 −0.32 ± 1.86 −0.31 ± 2.20 mth min [ • C] −0.02 ± 2.04 −0.23 ± 1.73 +0.11 ± 2.15 lakes with a delay of ∼ 1 month in the 1 • C warming day is much greater for deep high latitude lakes, then for low latitude or shallow lakes.Great Slave (62 • N and 41 m in depth) and Great Bear (66 • N and 72 m in depth) show a 28 day and 32 day delay in 1 • C warming day and a JAS LSWT decrease of 4.26 and 3.40 • C, while delays of 29 days and 32 days in the 1 • C warming day for Winnebago (44 • N) and Khanka (45 • N) both with depths of 5 m, showed only a small JAS LSWT decrease of ∼ 0.1 • C. In Fig. 17, the lake-mean depth of the 21 trial lakes is plotted against latitude.The relationship between the depth and latitude of the lakes and the decrease in the JAS LSWT (presented as the decrease in the JAS LSWT, per week of later 1 • C warming day, • C week −1 ) is shown in this figure, by use of coloured circles.This figure shows that for deep high latitude lakes the decrease in the JAS LSWT, is more pronounced than for shallow low latitude lakes.

Figure 16 .
Figure16.Lake surface water temperatures (LSWTs) for Great Bear (66 • N, 121 • W) and Great Slave (62 • N 114 • W) modelled with low snow and ice albedo (default albedo, α1: snow and white ice = 0.60 and melting snow and blue ice = 0.10) and high albedo (α2: snow and white ice = 0.80 and melting snow and blue ice = 0.60) demonstrating that the higher snow and ice albedo delays the 1 • C warming day, causing a lower July-August-September LSWT.

Figure 17 .
Figure17.The relationship between latitude and lake-mean depth of the 21 trial seasonally ice-covered lakes and the decrease in the July-August-September (JAS) lake surface water temperature (LSWT) caused by the later 1 • C warming day (as a result of using a high albedo, α2: snow and white ice = 0.80 and melting snow and blue ice = 0.60 in place of the default albedo, α1: snow and white ice = 0.60 and melting snow and blue ice = 0.10).The changes in the JAS LSWT, presented as the decrease in the JAS LSWT, per week of later 1 • C warming day, • C week −1 , are categorised by coloured circles.This figure indicates that high latitude and deep lakes show a larger decrease in the JAS LSWT per week of later 1 • C warming day, signifying that the LSWTs of these lakes are more responsive to changes in the 1 • C warming day, than low latitude and shallow lakes.

Figure 18 .
Figure 18.Comparison of lake-bottom temperatures during the stratification period, obtained from FLake model run using perpetual hydrological year, 2005/06 (Kirillin et al., 2011) and the monthly minimum climatology lake surface water temperature (LSWT) observations from ARC-Lake, for 14 deep (> 25 m) nonice-covered lakes (55 • S to 40 • N).The monthly minimum observed LSWTs have a ∼ 1 : 1 relationship with the lake-bottom temperatures during the stratification period.

Figure 19 .
Figure19.The lake mean depth vs. the modelled effective depth for 244 tuned lakes.Colour coding illustrates the effective depth factors.The average lake depth for each effective depth factor used in the tuning process is also given (insert).This figure demonstrates that deeper lakes are tuned to a shallower effective depth and shallower lakes to a deeper effective depth.

Table 3 .
Relationship between the lake surface water temperature (LSWT)-regulating properties and metrics, showing the equations for determining the daily mean absolute difference (MAD) and the July-August-September (JAS) LSWT mean difference.Units are given in square brackets.

Table 6 .
• N/S and no scaling (u1) to non-ice-covered lakes < 35 • N/S.Units are given in square brackets.Comparison of metric results for seasonally ice-covered lakes: 135 lakes tuned using the initial tuned set-up for seasonally icecovered lakes (Table2), 25 lakes tuned with the modified set-up (Table2), all lakes, and trial lakes.The spread of differences across lakes is defined as 2σ .The metrics are explained in Table4.Units are given in square brackets.

Table 7 .
Comparison of tuned model results for saline, freshwater, high and low altitude seasonally ice-covered lakes, with the spread of differences across lakes, 2σ .The metrics are explained in Table4.Units are given in square brackets.

Table 8 .
Comparison of tuned metric results for saline, freshwater and high and low altitude non-ice-covered lakes, with the spread of differences across lakes, 2σ .The metrics are explained in Table4.Units are given in square brackets.

Table 9 .
The fraction of observed inter-annual variance detected in the model.Maximum and minimum LSWT is used for non-ice-covered lakes (inter max and inter min ), while July-August-September (JAS) LSWTs are used for seasonally ice-covered lakes, (inter jas ).This table highlights that where the observed inter-annual variance is low, the proportion of variance detected in the model is also low (high altitude seasonally ice-covered lakes and tropical lakes).Units are given in square brackets.
lakes with greater observed variance have a greater portion of the variance detected in the model.For high altitude seasonally ice-covered temperate lakes, the fraction of the observed JAS LSWT inter-annual variance explained by the tuned model is considerably less (inter jas = 0.21) than for range of monthly LSWTs for non-ice-covered lakes, explain 0.38 and 0.36 (p < 0.0005) of the variation in var max and var min , with lakes of a low annual range (high altitude and tropical lakes), showing less inter-annual variance.This supwww.geosci-model-dev.net/9/2167/2016/Geosci.Model Dev., 9, 2167-2189, 2016 • C, has an earlier 1 • C cooling day, by 7.6 days.The deepest lake in the trials, Lake Hovsgol, average depth of www.geosci-model-dev.net/9/2167/2016/Geosci.Model Dev., 9, 2167-2189, 2016