Assessing degradation and recovery pathways in lakes impacted by eutrophication using the sediment record

Bennion, Helen; Simpson, Gavin Leslie; Goldsmith, Ben J.

doi:10.3389/fevo.2015.00094

ORIGINAL RESEARCH article

Front. Ecol. Evol., 11 August 2015

Sec. Paleoecology

Volume 3 - 2015 | https://doi.org/10.3389/fevo.2015.00094

This article is part of the Research Topic Using paleolimnology for management and restoration of lakes View all 12 articles

Assessing degradation and recovery pathways in lakes impacted by eutrophication using the sediment record

$\r\nHelen Bennion*$ Helen Bennion¹^*

Gavin L. Simpson^2,3

Ben J. Goldsmith¹

¹Department of Geography, Environmental Change Research Centre, University College London, London, UK
²Institute of Environmental Change and Society, Faculty of Science, University of Regina, Regina, SK, Canada
³Department of Biology, University of Regina, Regina, SK, Canada

Efforts to restore enriched lakes have increased yet there remains uncertainty about whether restoration targets can be achieved and over what timescale. Paleoecological techniques, principally diatom analyses, were used to examine the degree of impact and recovery in 12 European lakes subject to eutrophication and subsequent reduction in nutrient loading. Dissimilarity scores showed that all sites experienced progressive deviation from the reference sample (core bottom) prior to nutrient reduction, and principal curves indicated gradual compositional change with enrichment. When additive models were applied to the latter, the changes were statistically significant in 9 of the 12 sites. Shifts in diatom composition following reduction in nutrient loading were more equivocal, with a reversal toward the reference flora seen only in four of the deep lakes and one of the shallow lakes. Of these, only two were significant (Lake Bled and Mjøsa). Alternative nutrient sources seem to explain the lack of apparent recovery in the other deep lakes. In three shallow lakes diatom assemblages were replaced by a community associated with lower productivity but not the one seen prior to enrichment. Internal loading and top down control may influence recovery in shallow lakes and climate change may have confounded recovery in several of the study sites. Hence, ecosystem recovery is not simply a reversal of the degradation pathway and may take several decades to complete or, for some lakes, may not take place at all. By assessing ecological change over a decadal to centennial timescale, the study highlights the important role that paleolimnology can play in establishing a benchmark against which managers can evaluate the degree to which their restoration efforts are successful.

Introduction

Most lakes throughout the world have been modified to some extent by human activity. Excessive nutrient and organic matter loading resulting in eutrophication has affected numerous waterbodies, most notably since the mid-twentieth century associated with the intensification of agriculture and expansion of populations connected to sewage treatment works (Joye et al., 2006; Battarbee et al., 2011). The consequent high algal biomass leads to filtration problems for the water industry, oxygen depletion, recreational impairment, loss of biodiversity, fish mortality, and decline or loss of submerged plants (Smith et al., 1999).

Efforts to restore enriched systems have increased over the last few decades and there are now numerous examples of lakes in recovery (Jeppesen et al., 2005; Verdonschot et al., 2013). Point-source control at sewage treatment works has been particularly effective at reducing external nutrient loads but nutrient pollution from diffuse agricultural sources has proved more difficult to control as it is dispersed over large areas (Carpenter et al., 1998; Schoumans et al., 2014). Nevertheless, restoration schemes that promote use of buffer strips, good agricultural practice, and wetland regeneration have all contributed to the reduction of nutrient loading from agricultural sources (Sharpley et al., 2000). In deep, well flushed lakes, eutrophication is often reversed by the reduction in phosphorus (P) inputs alone, such as in Lake Washington, USA, where P concentrations fell dramatically, phytoplankton biomass declined and there were sustained increases in transparency following effluent diversion and treatment (Edmondson and Lehman, 1981). However, in shallow lakes internal P loading can delay recovery and external P reduction is often combined with other management measures such as dredging or biomanipulation (Søndergaard et al., 2007; Jeppesen et al., 2012), and increasingly geo-engineering techniques which use P-capping agents (Spears et al., 2013; Zamparas and Zacharias, 2014).

Whilst there are many individual success stories, there remains considerable uncertainty about whether restoration targets can be achieved and over what timescales one might expect to see improvement. Recovery may be a slow process as biotic communities tend to exhibit hysteresis and time-lags, and thus ecosystems take time to re-adjust to reduced stress (e.g., Yan et al., 2003; Johnson and Angeler, 2010). In an analysis of long-term datasets from 35 restored lakes, Jeppesen et al. (2005) showed that internal nutrient loading delayed recovery, but in most lakes a new equilibrium for total P (TP) was reached after 10–15 years. Furthermore, new pressures, especially from global warming, may counter restoration strategies. Climate change in combination with land use changes is anticipated to cause increased nutrient loading in lakes, and may increase the frequency and intensity of harmful algal blooms (Jeppesen et al., 2010, 2014). Longer growth seasons, higher water temperature, and more turbid conditions are likely to amplify eutrophication problems (Jeppesen et al., 2010; Moss et al., 2011). Thus the expectation that ecosystems can be returned, following remediation efforts, to conditions prior to enrichment may be a naive one and managers and policy makers may have to accept that “shifting baselines” will limit the ability to meet restoration targets (Duarte et al., 2009; Bennion et al., 2011a; Battarbee et al., 2012). Indeed the concept of “novel ecosystems” has been introduced to describe ecosystems in which the species composition and/or function have been completely transformed from the historic system (Hobbs et al., 2009). It has been argued that in these systems a refocusing of conservation and restoration practices away from existing or historical assemblages may be required though most ecologists advocate management based on the same values that we want from historical ecosystems (Higgs et al., 2014; Corlett, 2015).

Legislative programmes are now in force to reduce pollution and restore aquatic ecosystems to good health in many regions of the world. In Europe, the Water Framework Directive (WFD) with its aim to restore waters to at least good status, has increased the need for effective restoration programmes for all lakes (European Union, 2000). Within the WFD, ecological status is based on the degree to which present day conditions deviate from those expected in the absence of significant anthropogenic influence, termed reference conditions. Consequently there has been a wave of research aimed at defining reference conditions and development of tools for estimating deviation from them. Lake sediment analysis provides unique insights into the history of freshwater ecosystems giving evidence for the nature and timing of ecosystem change, and providing a record of human impact that can be indispensible in developing strategies for ecosystem management (Bennion et al., 2011a). Palaeoecological methods can reveal pre-impact conditions and identify any signs of recovery, and have played a key role in the WFD (Bennion and Battarbee, 2007), particularly in determining pre-enrichment reference conditions and degree of eutrophication (Bennion et al., 2004). Diatom records have proved especially valuable in this respect, largely due to their sensitivity to shifts in trophic status (Bennion and Simpson, 2011; Bennion et al., 2011b). As many restoration programmes progress, there is great potential to employ a combination of limnological and sediment records to track recovery using the pre-eutrophication baseline as a benchmark (Battarbee et al., 2005).

This paper employs palaeoecological techniques to examine the degree of impact and recovery in 12 European lakes that have been subject to eutrophication. Changes in the diatom assemblages, both community composition and diatom-inferred TP (DI-TP) concentrations, in sediment cores from the study lakes are assessed in response to enrichment and subsequent reduction in nutrient loading. The extent to which the diatom assemblages are approaching or deviating from reference conditions is explored using ordination and dissimilarity scores and the identification of the floristic changes is assessed statistically using additive models. The key questions being addressed are: (i) Do the observed changes reflect degradation and recovery?, (ii) Is the recovery pathway simply a reversal of the degradation pathway?, and (iii) Can the lake sediment record be used to track degradation and recovery and thus inform management?

Materials and Methods

Study Sites

The 12 study sites are located in European lowland catchments and represent a range of lake types in terms of lake area, depth, and trophic status (Table 1). The lakes are within the temperate climatic zone, located in six countries across Europe from Slovenia to Norway, spanning a latitudinal range from 46.4 to 60.8°N. Most of the lakes lie in lowland (< 200 m) catchments, with the exception of Lake Bled in Slovenia which lies at higher elevation. The dataset covers a wide range of lake surface areas from < 1 to 362 km² though is comprised mostly of small lakes. Similarly, lakes span a broad range in maximum water depth from 1.5 to 453 m. For data analyses, the sites have been classed as either deep, stratifying (seven lakes) or shallow, non-stratifying (five lakes), in order to explore whether these lake types respond differently to nutrient reduction measures. While Loch Leven has a deep basin of 25.5 m, the majority of the lake is shallow with a mean depth of 4 m, and as the core was collected in a water depth of 4 m, this lake is classed as shallow, non-stratifying.

TABLE 1

Table 1. Summary characteristics of the 12 study sites.

The lakes cover a wide range of current (annual mean) TP concentrations from 4 to 476 μg L⁻¹, spanning the full trophic gradient from oligotrophic to hypertrophic conditions (Table 1). The catchments are largely productive with nutrient loading from either point sources such as sewage treatment works and/or diffuse sources from agriculture. All of the sites have experienced eutrophication within the last ~100 years and have either seen a reduction in external nutrient loading from sewage treatment works since the 1970s and 1980s or have received less diffuse pollution as a result of changes in land use and management (Table 1).

Sediment Core Collection and Analyses

A sediment core was collected from the open water area of each lake as part of several different previous studies and, therefore, coring methods and analytical resolution vary from site to site. At least 10 samples from each core spanning the last ~200 years were analyzed for diatoms and, while there were typically >5 samples representing the post-restoration period, in a few cases there were only three samples available. Diatom analysis was carried out using standard methods (Battarbee et al., 2001). A minimum of 300 valves were counted from each sample using a research quality microscope with a 100x oil immersion objective and phase contrast. Krammer and Lange-Bertalot (1986–1991) was the principal flora used in identification. The diatom data are expressed as percentage relative abundances.

Chronologies for the cores included in this study were determined using radiometric methods. Selected sediment samples were analyzed for ²¹⁰Pb, ²²⁶Ra, ¹³⁷Cs, and ²⁴¹Am by direct gamma assay using the methods of Appleby et al. (1986). ²¹⁰Pb chronologies were calculated using either Constant Rate of Supply (CRS) or Constant Initial Concentration (CIC) dating models (Appleby and Oldfield, 1978), based on the method best suited to the data.

Data Analysis

The degree of floristic change in the diatom assemblages between the bottom sample and every other sample in the core was assessed using the squared chord distance (SCD) dissimilarity coefficient (Overpeck et al., 1985) implemented in C2 (Juggins, 2003). The dissimilarity scores range from 0 to 2 whereby 0 indicates that two samples are exactly the same and 2 that they are completely different. This provides a measure of deviation from the reference assemblage. For seven cores, the records extend back to 1800 AD and, therefore, represent a time period prior to major industrialization and/or agricultural intensification (Bennion and Simpson, 2011; Battarbee et al., 2011). For the remaining sites, the cores do not extend back this far but do cover the pre-enrichment period. It is generally considered that a shift from benthic to planktonic dominance occurs with eutrophication (Vadeboncoeur et al., 2003) and such a shift has been noted in the diatom assemblages in several palaeoecological studies (e.g., Battarbee, 1978; Sayer et al., 1999). Hence, the percentage of planktonic taxa vs. non-planktonic taxa was calculated for each sample to assess whether this provides a useful metric for assessing diatom response to degradation and recovery. A simple classification was used such that all centric diatoms, Asterionella formosa, Tabellaria flocculosa, and Fragilaria crotonensis were classed as planktonic and all others as non-planktonic taxa. This was considered sufficient for revealing the broad habitat shifts.

The first axis scores of a principal components analysis (PCA) or a correspondence analysis (CA or detrended CA - DCA) often fail to capture a long or dominant gradient such as we might expect to be present in temporally-ordered data with progressive change in abundance or composition of organisms. As a result, the time series of axis one scores may be a poor summary of compositional change (Simpson and Birks, 2012). Principal curves (Hastie and Stuetzle, 1989; De'ath, 1999; Simpson and Birks, 2012; Simpson, unpublished; PrC) is a non- or semi-parametric alternative to PCA, CA, DCA, etc. that is particularly suited to the identification of single or dominant gradients within a sediment core sequence. Here we use PrCs to summarize the timing of the major compositional changes in the diatom profiles for the 12 study lakes. We follow Simpson and Birks (2012) and allow the spline degrees of freedom to vary between species. PCA axis 1 was used as the starting curve in the principal curve fit for all sites except Mill Loch and Mjøsa, where a CA axis 1, and Gjersjøen, where PCA axis 2, were used respectively. These changes were needed to achieve a satisfactory fit with simple species response curves along the fitted PrC.

The PrC scores were extracted for each fit and arranged in time order. To determine if statistically significant change in composition (e.g., eutrophication and subsequent recovery) could be identified, we modeled the time series of PrC scores using additive models, with a continuous time first-order autoregressive (CAR(1)) process for the residuals to account for the lack of independence between observations (Simpson and Anderson, 2009). Note that here we use sample age as the sole covariate in the model and therefore, unlike the examples in Simpson and Anderson (2009), we simply wish to estimate the potentially non-linear trend in the PrC scores for each lake, to avoid over-interpretation of the time series of scores. Many of the cores were sampled with strongly-varying density in time, and often considerably more samples were available in the most recent sediments than the reference period. Therefore, we chose to place knots at the deciles of the distribution of sample ages for each core. This allowed the trend splines to adapt to the data in the regions of the cores where more observations were available. Fitting the additive model to the Kiełpińskie data required additional steps beyond those outlined here and in Simpson and Anderson (2009). It was not possible for the optimisation process to automatically distinguish between a wiggly trend or strong autocorrelation and the resulting model fit was highly uncertain owing to poor estimation of both the trend and the CAR(1) term in the covariance matrix of the residuals. To achieve an acceptable fit to the Kiełpińskie data we fixed the CAR(1) parameter, φ, equal to 0.5, allowing for moderate to strong residual autocorrelation, and limited the dimension of the spline basis functions to k = 8 (k = 10 in the other fits). This achieved a balance between accounting for residual autocorrelation whilst not over-fitting the data.

Approximate significance of the fitted trends can be achieved through the usual methods of statistical inference for additive models (Wood, 2006), however, this only provides a test against the simple null hypothesis of no change in PrC scores with time. Here, our interest is also on where in time the additive models suggest that compositional change takes place if the null hypothesis is rejected. In a linear regression, the slope of the regression line is the first derivative of that curve and, given the standard error of the estimate of the slope, one can determine if the slope differs from 0 (i.e., no effect). With the additive model, the slope of the fitted trend is potentially changing continuously over the time series of PrC scores, and as a result we do not have a single measure of departure from zero slope. Instead, we use the first derivative of the fitted trend spline evaluated at a set of regularly spaced time points over the interval covered by each sediment core.

To estimate the first derivative of the fitted trend spline we use the method of finite differences, in which we predict from the model at the set of regularly-spaced time points and at a second set, shifted relative to the first by a very small amount. The rate of change in the predictions between the original and shifted points is an approximation to first derivative. How well the finite difference method approximates the unknown first derivative of the spline is governed by size of the shift, with smaller values producing more accurate approximations; we used a value of 0.000001 (of a year) as the shift but the results are not sensitive to a range of values around the chosen value. With each first difference estimate we also obtained its standard error and then formed a pointwise 95% confidence interval on estimate. Where this confidence interval excludes 0 (no change, zero slope) we conclude that significant change in PrC scores is observed for that time point. We indicate these periods of statistically significant change using thicker sections of the fitted trend when plotted. Color is also used to convey meaning; blue indicates significant decrease in PrC scores whilst red indicates significant increase. In some plots, we have negated the PrC scores such that enrichment in the diatoms is associated with smaller PrC scores, and less enriched periods with larger PrC scores. This is justified as the PrC score is defined as the arc length along the fitted principal curve from one arbitrarily-chosen end of the curve (Hastie and Stuetzle, 1989). Hence, the sign on the scores is arbitrary and negating the scores for some lakes improves comparison of the extracted PrC scores and the additive models fitted to them.

A diatom-TP transfer function was applied to the diatom data to reconstruct the trophic status of each site using either: (1) a training set of 56 relatively large, deep lakes (>10 m maximum depth) from Scotland, Northern Ireland, Cumbria, southern Norway, and central Europe with annual mean TP concentrations ranging from 1 to 73 μg TP L⁻¹ and a median value of 22 μg TP L⁻¹ (Bennion et al., 2004); the best model was generated with simple weighted averaging and inverse deshrinking (ter Braak and van Dam, 1989); this model has a coefficient of determination (r²) between observed and inferred values of 0.75 and a root mean squared error of prediction (RMSEP based on the jack-knifing cross validation method) of 0.25 log₁₀ μg TP L⁻¹, or 2) a Northwest European training set of 152 relatively small, shallow lakes (< 10 m maximum depth) with annual mean TP concentrations ranging from 5 to 1190 μg TP L⁻¹ and a median value of 104 μg TP L⁻¹, and an r² of 0.91 and a RMSEP of 0.21 log₁₀ μg TP L⁻¹ for the weighted averaging partial least squares two-component (WA-PLS2) model (Bennion et al., 1996). The most appropriate transfer function was used in each case, i.e., the training set containing the greatest percentage of taxa present in the fossil samples. The DI-TP values for the surface sediment samples were generally in good agreement with the current measured TP concentrations of the study lakes. All reconstructions were implemented using C2 (Juggins, 2003).

The samples of each core were projected passively into a PCA of samples from the modern diatom-TP training sets described above, forming so called timetrack plots. These plots allow the direction of floristic change at each site to be visualized. The sample and species scores were plotted in ordination biplots to illustrate the degree to which the recovery trajectories follow back along the enrichment trajectories. For deep lakes only species that are present in at least two sites and with maximum abundance greater than or equal to 2% are shown and, for shallow lakes, these values were five occurrences and 5% abundance, respectively. These choices were purely for aesthetic reasons to improve the clarity of the plot whilst retaining labeling for the most abundant or frequently occurring taxa. Arrows have been added to the plots to illustrate the enrichment and post-restoration trajectories.

Unless stated otherwise, all analyses were conducted using R (version 3.1.3, R Core Team, 2015) with several additional packages: vegan (version 2.2-1, Oksanen et al., 2014) was used to fit the ordinations, analog (version 0.16-0, Simpson, 2007; Simpson and Oksanen, 2014) was used to fit the principal curves and time tracks, and additive models were fitted using the mgcv package (version 1.8.4, Wood, 2004, 2006, 2011). Additional R functions written by GLS (based on suggestions by Simon Wood, pers. comm.) were used to evaluate the first derivative of the trend splines and form the pointwise confidence interval. R scripts implementing the analyses and reproducing the figures are available online from https://github.com/gavinsimpson/bennion-frontiers-2015.

Results

Dissimilarity Scores

The dissimilarity scores between core bottom and other samples in each core indicate that all sites have experienced deviation from reference condition (core bottom sample) over the period represented by the cores (Figure 1). All sites, with the exception of Kiełpińskie and Lidzbarskie where the patterns are less clear, exhibit progressive deviation from the reference sample during the period prior to nutrient reduction, indicating gradual compositional change with enrichment. The diatom assemblages of some sites, most notably the deep lakes, show signs of returning toward the reference flora following reduction of nutrient load. This is most apparent in Lake Bled and Mjøsa, and to a lesser extent in Gjersjøen. Nonetheless, most are still far from reference condition with high dissimilarity scores ranging from 0.38 to 1.57 between the core top and bottom samples.

FIGURE 1

Figure 1. Dissimilarity scores for the 12 study sites. The degree of floristic change in the diatom assemblages between the bottom sample and every other sample in the core; timing of first nutrient reduction is shown by the arrow.

Percentage Plankton

In four of the five shallow lakes (Barton Broad, Loch Leven, Llangorse Lake, and Marsworth Reservoir) the % plankton increases with enrichment but does not decline during the recovery phase (Figure 2). In the other shallow lake, Bosherston Central Lake, the % plankton stays low throughout the record. In the deep lakes % plankton was high throughout the cores (generally >60%) but in Esthwaite Water, Gjersjøen, Mill Loch and Mjøsa slight increases in the planktonic component were observed with enrichment. Only in Mjøsa, and to a lesser extent in Esthwaite Water, was a slight decline in % plankton seen in the recovery period.

FIGURE 2

Figure 2. Percentage plankton in the 12 study sites. Timing of first nutrient reduction is shown by the arrow.

Ordination and Transfer Functions

The PrC analysis shows that all sites experienced gradual yet unidirectional shifts in PrC scores during the eutrophication phase, occurring for example in Marsworth Reservoir from 1960 to 1980, in Lake Bled from 1940 to 1970 and in Gjersjøen from 1960 to 1980 (Figure 3). However, notwithstanding the slight offset in timing, the scores move in the reverse direction at only four sites (Marsworth Reservoir, Lake Bled, Gjersjøen, and Mjøsa) following reduction in nutrient loading (Figure 3). At the remaining sites, there is no clear pattern in the direction of the PrC scores following restoration.

FIGURE 3

Figure 3. Principal curve scores for the 12 study sites. Timing of first nutrient reduction is shown by the arrow.

The additive model fits to the PrC scores indicate change in one direction to some degree in all sites except Gjersjøen, where the rapid shifts seen across relatively few data points are not well captured by the model and hence a flat line (no change, with strong residual autocorrelation) is the best fitting model (Figure 4). The changes are significant (as indicated by the blue lines in Figure 4) in all but Marsworth Reservoir and Lidzbarskie. At Marsworth Reservoir, some of the patterns in the data are explained as temporal autocorrelation rather than as trend. At Lidzbarskie, there is considerable scatter in the points and, once temporal autocorrelation is allowed for, the uncertainty on the fitted trend is sufficiently large that a statistically significant pattern cannot be extracted from the data series. For the deep lakes, the first significant decrease signaling eutrophication occurs at c.1850 in Mjøsa and Mill Loch, c.1910 in Kiełpińskie, c.1930 in Esthwaite Water, and c.1940 in Lake Bled (Figure 4B). For the shallow lakes, the first significant decrease occurs at c.1850 in Loch Leven, c.1910 in Barton Broad, c.1940 in Llangorse Lake, and in the early 1980s in Bosherston Central Lake (Figure 4A). A statistically significant recovery trend is observed only in Lake Bled and Mjøsa (as indicated by the red line in Figure 4B), starting at around the early-1980s and mid-1970s, respectively, and in both cases is approximately coincident with restoration measures although there is a slight offset, with the shifts in the core pre-dating the introduction of measures by around 5 years. At Kiełpińskie, significant compositional change is evident from c.1945 to 1960 (red line in Figure 4B) but this reflects modest changes in Aulacoseira spp. and Cyclotella spp. and does not appear to be associated with recovery from eutrophication as Stephanodiscus parvus, typically found in nutrient-rich waters, continues to increase and become dominant during this period.

FIGURE 4

Figure 4. Principal curve fits with derivatives for (A) shallow lakes and (B) deep lakes. Blue and red lines indicate significant decrease and significant increase in PrC scores, respectively. Dashed lines are 95% confidence intervals. Timing of first nutrient reduction is shown by the arrow.

The diatom transfer functions infer an increase in TP concentrations in seven of the study lakes during the enrichment period (Figure 5). Of the shallow lakes, a clear signal was not seen in Barton Broad, Bosherston Central Lake or Llangorse Lake, and of the deep lakes a clear increase was not apparent in Gjersjøen or Kiełpińskie. However, the diatom transfer functions infer a decline in TP concentrations following a reduction in nutrient loading at 11 of the 12 study lakes, the exception being Kiełpińskie where a decrease in DI-TP concentrations is not clearly seen (Figure 5). This suggests that at these 11 sites there have been compositional changes toward taxa associated with lower nutrient concentrations following the nutrient reductions. In the case of Kiełpińskie, the shifts in the diatom assemblages were subtle and have resulted in no major change in DI-TP values in recent years.

FIGURE 5

Figure 5. Diatom-inferred TP (DI-TP) reconstructions for the 12 study sites. Timing of first nutrient reduction is shown by the arrow.

The PCA of the deep lakes training set revealed that axis 1 (PC1: eigenvalue 10.43) and axis 2 (PC2: eigenvalue 7.76) explained 16% and 12% of the variance in the species data, respectively. These values are relatively large for diatom training sets especially those with very many taxa as here. When the deep lake cores are plotted passively on a PCA of the large, deep lakes training set samples (Figure 6) the core samples generally move from the lower right of the plot toward the upper left during the enrichment period. A floristic reversal is most clearly seen in Lake Bled (A), Gjersjøen (C) and Mjøsa (G), where samples move back toward the lower right following a reduction in nutrient loading, and to a lesser extent in Mill Loch (F). However, the additive models indicate that this reversal is statistically significant only in the former two sites. This reverse pattern is not seen in Esthwaite Water (B). Nor is a clear pattern observed for the two Polish lakes: at Kiełpińskie (D) samples move from right to left but there is no subsequent reversal, and at Lidzbarskie (E) there is no clear direction of change. The core trajectories reflect changes in the composition of the diatom flora with taxa associated with lower nutrient concentrations located on the right of the diagram (e.g., Achnanthes spp., Brachysira spp., Cymbella spp., Eunotia spp., oligotrophic Cyclotella spp., and Tabellaria flocculosa), those more typically found in waters with intermediate nutrient concentrations located on the upper left (e.g., Aulacoseira subarctica, Asterionella formosa, Fragilaria crotonensis, and Cyclotella radiosa) and those commonly seen in nutrient-rich conditions located on the lower left of the plot (e.g., Aulacoseira granulata, Cyclostephanos dubius, Stephanodiscus hantzschii, Stephanodiscus parvus) (Figure 7). The samples for Kiełpińskie and Lidzbarskie are compressed into a small space (Figure 6) owing to the dominance of Stephanodiscus spp. and Cyclostephanos spp. throughout the records.

FIGURE 6

Figure 6. Deep lake cores plotted passively on a PCA of training set samples. The direction of change over time is shown by the arrows. (A) Lake Bled, (B) Esthwaite Water, (C) Gjersjøen, (D) Kielpinskie, (E) Lidzbarskie, (F) Mill Loch, (G) Mjøsa.

FIGURE 7

Figure 7. Species plot on a PCA of training set samples for the deep lake cores. (A) All taxa present in at least two sites and with maximum abundance greater than or equal to 2%, (B) zoomed in section as indicated by the inset in (A). See Supplementary Material, Appendix 1 for diatom codes and names.

The PCA of the shallow lakes training set revealed that PC1 (eigenvalue 9.93) and PC2 (eigenvalue 7.15) explained 15% and 11% of the variance in the species data, respectively. When the shallow lake cores are plotted passively on a PCA of the shallow lake training set samples (Figure 8) the core samples of all lakes, with the exception of Bosherston Central Lake (B), move from the left of the plot toward the right during the enrichment period. A clear floristic reversal is apparent only at Marsworth Reservoir (E) following nutrient reduction, while a slight move back toward the left of the diagram is seen at Loch Leven (C), though additive models indicate that the reversal is not statistically significant in either case. At Barton Broad (A), Bosherston Central Lake (B), and Llangorse Lake (D), the upper core samples move to a new position within the ordination space but do not obviously track back along the enrichment trajectory. As for the deep lakes, the core sample shifts during the eutrophication phase largely reflect a move from taxa associated with relatively nutrient poor conditions located to the left of the diagram (e.g., Achnanthes spp., Brachysira spp., Cymbella spp., oligotrophic Cyclotella spp., and Tabellaria flocculosa) to those taxa typically found in nutrient rich waters located to the right of the plot (e.g., Cyclotella meneghiniana, Cyclostephanos dubius, Cyclostephanos tholiformis, Stephanodiscus hantzschii, Stephanodiscus parvus) (Figure 9). The benthic Fragilaria spp. which are often abundant in shallow lakes, were positioned in the upper left of the diagram.

FIGURE 8

Figure 8. Shallow lake cores plotted passively on a PCA of training set samples. The direction of change over time is shown by the arrows. (A) Barton Broad, (B) Bosherston Central Lake, (C) Loch Leven, (D) Llangorse Lake, (E) Marsworth Reservoir.

FIGURE 9

Figure 9. Species plot on a PCA of training set samples for the shallow lake cores. (A) All taxa present in at least five sites and with maximum abundance greater than or equal to 5%, (B) zoomed in section as indicated by the inset in (A). See Supplementary Material, Appendix 1 for diatom codes and names.

Discussion

Diatom Response to Changes in Nutrient Loading

One of the challenges for ecologists wishing to track environmental change is to find biological indicators that are sufficiently sensitive to the pressure gradient of interest. Here, a range of diatom metrics were explored as diatoms are sensitive to changes in water quality and are particularly good indicators of lake nutrient concentrations (Hall and Smol, 2010). The most striking changes were observed in diatom composition and were effectively summarized by the PrC and dissimilarity (SCD) scores. The data demonstrate that progressive deviation from the reference condition (here defined as the assemblage at the bottom of the core) occurred at all sites during the eutrophication phase although the trends were statistically significant at only nine of the 12 sites. The diatom shifts were gradual rather than abrupt, reflecting a process of relative decline in taxa associated with low nutrient concentrations and their replacement with taxa typically found in more nutrient-rich waters. The ordination plots illustrate that whilst the reference conditions of the 12 study lakes are site specific there are some common patterns in compositional change with shifts from a flora composed of Achnanthes spp., Brachysira spp., Cymbella spp., Eunotia spp., oligotrophic Cyclotella spp., and Tabellaria flocculosa to one composed of Aulacoseira subarctica, Aulacoseira granulata, Asterionella formosa, Fragilaria crotonensis, and Cyclotella radiosa as enrichment progresses, and in the most nutrient-rich cases, to an assemblage composed of small centric taxa such as Cyclostephanos dubius, Stephanodiscus hantzschii, and Stephanodiscus parvus. These same shifts have been observed in numerous European lakes during periods of increased nutrient loading (e.g., Anderson, 1997; Lotter, 1998, 2001; Bennion et al., 2004, 2011b) and, therefore, provide a useful indication of ecological change associated with eutrophication.

The shifts in diatom composition following reduction in nutrient loading are more equivocal. A clear reversal toward the reference flora is seen only in three of the deep lakes, Lake Bled, Gjersjøen, and Mjøsa, and to a lesser extent in Mill Loch, a deep lake, and Marsworth Reservoir, a shallow lake, in terms of both the dissimilarity and direction of PrC scores. However, the model fits indicate that the recovery trend is only statistically significant at Lake Bled and Mjøsa. It is noteworthy that the three lakes that show the clearest reversal toward the reference flora are those with the lowest mean TP concentrations and are the deepest lakes in the dataset, supporting observations elsewhere that in deep lakes eutrophication is often reversed by the reduction in P inputs alone (Jeppesen et al., 2005). Nevertheless, as for the degradation phase, the compositional changes after nutrient reduction are gradual rather than sudden suggesting that ecological recovery may take several years to decades to complete. Indeed even these five lakes do not exhibit a return to the pre-enrichment flora over the 20–30 year period since remedial measures were introduced, as dissimilarity scores between the core bottoms and tops remain relatively high, and the core trajectories illustrate that the recent assemblages have not yet returned to those observed in the lower cores. The data for these five lakes suggest that, whilst the diatoms have responded to nutrient reduction and are heading back along the eutrophication pathway, they still have some way to go before they reach reference condition.

For the remaining lakes, the diatom response during the “recovery” period is more difficult to discern. At Esthwaite Water, P loading from a local sewage treatment works was reduced in 1986 but nutrients derived from catchment runoff, a fish farm established in 1981 and sediment P release have negated any potential reduction in lake nutrient concentrations (Bennion et al., 2000; Dong et al., 2012) and thus biological recovery is not apparent. For the two Polish lakes (Kiełpińskie, Lidzbarskie) which are relatively deep, stratifying waterbodies, the reversal toward former assemblages is less striking than for the other deep sites. This may be because there have been no specific restoration measures taken to reduce point sources of nutrients, and any decrease in nutrient loading is due to the reduction of fertilizer use and changes in land use in the catchments caused by the significant economic changes in the country in the early 1990s. Diatom response is, therefore, very recent and appears to be relatively subtle. The sediment records of these two lakes are also rather short, particularly so for Lidzbarskie where the core represents only 50 years. Analysis of a longer sediment core may be necessary to determine the pre-enrichment diatom assemblages and to provide a more comprehensive assessment of timing and degree of change at this site.

The Loch Leven data point to partial recovery of the diatom flora, namely an increase in Aulacoseira subarctica relative to Stephanodiscus taxa since the mid-1980s, as a result of a catchment management plan introduced in 1985 (Bennion et al., 2012). Long term datasets for the loch show that P concentrations have declined markedly but the trend was non-linear with a slight increase in the early 1990s caused by P recycling from the sediments (Carvalho et al., 2012). Additionally agriculture in the Leven catchment remains a significant diffuse source of nutrients to the loch as much of the land is used for arable farming, and rural septic tanks also contribute to the P load (May et al., 2012). The sediment record suggests that those diatom taxa lost during enrichment have not yet returned, most likely because nutrient concentrations remain too high (Bennion et al., 2012). Likewise, at Barton Broad there is little evidence of any recovery in the diatom assemblages. In spite of a substantial reduction in the amount of P entering the rivers from sewage treatment discharges and a consequent progressive decline in lake TP and chlorophyll a concentrations since the late 1970s and early 1980s, the reduction in epilimnetic TP was slow due to the continued release of P from the sediments (Phillips et al., 1999, 2005). Barton Broad, therefore, remains dominated by phytoplankton with almost no submerged macrophyte growth and it is perhaps not surprising that the diatoms show minimal response.

At the other two shallow lakes (Bosherston Central Lake and Llangorse Lake) there are compositional changes in the diatoms following remediation efforts but the diatom floras do not appear to revert back toward those seen prior to enrichment but rather move toward a different assemblage. At Llangorse Lake, the effluent from the local sewage treatment works was diverted from the lake in 1981, with a second smaller input diverted in 1992. While the lake appears to be recovering following the remedial measures taken (Bennion and Appleby, 1999), the diatom community in the upper samples is dominated by planktonic forms that were not previously abundant in the record (namely Aulacoseira subarctica, Aulacoseira ambigua, and Cyclotella radiosa) and has not yet returned to the Fragilaria spp. dominated assemblage seen prior to enrichment. One explanation for this is that algal productivity is N-limited rather than P-limited for most of the summer and internal P loading from the sediments in summer remains high (May et al., 2010). Furthermore, sub-surface flow has been identified at Llangorse Lake with potential to bring in nutrients from outside the immediate catchment (May et al., 2010). At Bosherston Central Lake monitored total phosphate concentrations have exhibited a decline from the early 1980s following various interventions including diversion of sewage since 1984 and construction of a bypass pipeline in 1992 (Davidson et al., 2002). The principal change in the diatom assemblages from this time, most notably since the early 1990s, has been a shift in life-forms with an overall increase in epiphytic taxa relative to benthic taxa. Interpretation of the diatom species shifts is difficult owing to the subtle nature of the changes and the uncertainties regarding the factors which determine the composition of non-planktonic communities. Nevertheless, changes in the nutrient concentrations, shifts in habitat availability, and plant community structure and alterations in grazing pressures are all possible explanations (Davidson et al., 2002). The importance of both top down and bottom up mechanisms, the role of the trophic cascade, interactions between the pelagic and littoral environments, the potential for alternative stable states (Scheffer et al., 1993), and the numerous pathways of P recycling must all be considered in order to fully understand how shallow lake ecosystems might respond to changes in nutrient loading. It is important to bear this complexity in mind when attempting to interpret the findings for the shallow lakes in the present study.

The data suggest that the percentage of planktonic taxa may be a useful metric for tracking enrichment in shallow lakes as a shift toward higher percentage of plankton was apparent with increased nutrient loading in four of the five shallow sites. An increase was less obviously seen in the deep lakes where percentage of planktonic taxa was generally high throughout the records. The shift from benthic to planktonic production associated with eutrophication in shallow waterbodies has been well documented (e.g., Vadeboncoeur et al., 2003). Benthic algae often become light limited as planktonic forms become more abundant, and as submerged macrophytes are lost so too are potential habitats to support epiphytic taxa. However, a subsequent decline in the planktonic component of the diatom assemblages following nutrient reduction was not evident, indicating that there is a degree of “unhelpful resilience” (Standish et al., 2014) and hysteresis in the systems and the diatom flora does not automatically revert back to that seen prior to enrichment.

Increases in DI-TP were observed in seven lakes during the enrichment period while a decline in DI-TP was seen in 11 lakes following remediation. Hence this metric, which essentially reflects shifts in diatom composition, appears to have some potential for tracking recovery. However, several studies have highlighted the shortcomings of the transfer function technique in certain situations, and this is particularly well documented for shallow lakes where non-planktonic taxa dominate the diatom assemblages. Problems include the influence of factors such as light, substrate, and top-down factors in addition to water chemistry on the distribution of these taxa and their wide tolerance to nutrient concentrations, making them poor indicators of lake trophic status (e.g., Anderson et al., 1993; Bennion, 1995; Bennion et al., 2001; Sayer, 2001; Juggins et al., 2013). Nonetheless, in the absence of other techniques for hindcasting nutrient concentrations, inference models are likely to remain a valuable part of the lake manager's toolkit (Saulnier-Talbot, 2015).

Degradation vs. Recovery Pathways

Our palaeoecological data reveal that whilst in some cases the diatom recovery trajectories do appear to track back along the degradation pathway, in others and the shallow lakes in particular, either little sign of recovery is evident or the assemblages follow a new trajectory. It might be argued that a return to the pre-disturbance diatom flora is less likely than the establishment of a different but ecologically equivalent assemblage i.e., species with a similar autecology. This was seen in several of the shallow lakes (e.g., Bosherston Central Lake and Llangorse Lake) where there were compositional changes toward taxa associated with lower nutrient concentrations following nutrient reduction but these were not the same taxa as those present in the pre-enrichment samples.

Our data accord with the findings of other recovery studies in that ecosystem recovery is shown not simply to be a reversal of the degradation process. Duarte et al. (2009) examined four coastal systems demonstrating that they failed to return to the reference status upon nutrient reduction, offering alternative nutrient sources, internal loading, shifts in limiting nutrients, co-limitation effects of nutrients and light, and decreased filter-feeder activity as potential explanations for failure. Similarly, riverine communities do not necessarily show the anticipated and desired signs of improvement and recovery may lead to endpoints very different from the original undisturbed state (Palmer et al., 1997). Time lags associated with the release of legacy P (the surplus P stored in soils and sediments derived from past land use activities) is one plausible explanation for why P controls may not produce expected improvements in water quality (Jarvie et al., 2013; Sharpley et al., 2013). In UK lakes, Battarbee et al. (2014) demonstrated that recovery of diatom communities from acidification was limited when compared with the pre-acidification reference. Similarly, to the current study, in a few cases the floristic composition of recent samples was different from that observed during and before the acidification phase and, while the reasons for this are not yet clear, nutrient enrichment from atmospheric N deposition and/or climate change are potential factors confounding recovery (Battarbee et al., 2014). Lake ecosystems have been reported to follow convoluted trajectories following nutrient reduction, with internal loading, changes in food webs, the impacts of climate change, and 10–15-year time lags proposed as the causes for the complex lake trajectories observed (Jeppesen et al., 2005). In a palaeoecological study of six Swiss lakes, diatoms did not switch back to the ones characteristic of the early phase of enrichment despite reductions in nutrient loading (Lotter, 2001).

Most of the existing studies on recovery pathways are based on long-term datasets but for the majority of freshwater ecosystems monitoring activities are rather short-term and do not sufficiently account for long time periods required for restoration. The longer timeframe afforded by the sediment record thus lends itself well to studies of lake recovery and by extending back several decades or even centuries is valuable for defining the reference condition against which degree of recovery can be assessed. Nonetheless, palaeoecological data are not without their limitations and in this study it is perhaps the relatively low resolution of the data for the recovery period (i.e., at some sites only two or three samples correspond to the period since nutrient reduction) that is the greatest weakness. Studies are in progress whereby sites that were cored over 20 years ago are being repeat-cored and the sediments fine-sliced to track recovery over the last few decades at a higher resolution than was possible here. Error associated with the chronologies derived from radiometric dating, albeit small in recent sediments of typically only ±2–3 years, can cause difficulty in exactly pinpointing the timing of nutrient reduction in a core and may explain the slight offset between the management actions and diatom response observed in some of our records. Equally, discrepancy in the timings may arise due to changes in land use and management in the study catchments that are not documented and may have occurred prior to the more major and well documented management actions listed in Table 1. Furthermore, one has to allow for a good deal of approximation in the models and fits, accepting that the model will interpolate between adjacent time points and the strength of the fitted trend at any point is dependent upon the pattern in the data both prior and subsequent to any particular time point of interest. Caution must also be exercised when interpreting changes seen in surface sediment samples (upper 0–1 cm) in terms of recovery as the reduced period of decomposition experienced by recently deposited diatoms can result in these being over-represented (Sayer, 2001). Varved sequences, offering an annual resolution, lend themselves particularly well to studies of degradation and recovery (e.g., Chandler Rowell et al., 2015). A combination of long-term datasets and palaeolimnological approaches provides a particularly powerful tool for assessing timescales of ecological change (Battarbee et al., 2005; Bennion et al., 2012; Dong et al., 2012).

Management Implications and Factors Confounding Recovery

Several of the study lakes exhibit signs of ecological recovery in terms of reversal in their diatom assemblages but even in these cases the assemblages are still far from reference conditions as much as two to three decades since management measures were taken to reduce nutrient loads. Our findings are in accordance with other studies which report typical timescales of recovery for lakes of 10–20 years lakes (Jeppesen et al., 2005; Jones and Schmitz, 2009; Verdonschot et al., 2013). This has major implications for the WFD which requires waterbodies to be restored to at least good status, over the course of the next two river basin planning cycles in 2021 and 2027, in that the effects of any measures that have recently been introduced could take several decades to be seen. Perhaps even more importantly the data suggest that for some lake systems the assemblages following remedial action may not return back down the degradation pathway at all and, therefore, reference conditions are unlikely ever to be achieved. Factors such as climate change and diffuse, legacy, and emerging pollutants (Deblonde et al., 2011) that are difficult to control may limit the ability to meet restoration targets and managers may have to accept that species assemblages following remediation may not return to those of the past.

In most of our study lakes, the main point source of nutrients, principally P, has been the key focus of management efforts. However, in recent decades diffuse nutrient sources have become relatively more significant than urban wastewater pollution and losses from agricultural land are now the biggest challenge (Schoumans et al., 2014). There has been a growing literature on the need to reduce nitrogen (N) loads as well as P in order to reverse eutrophication (Galloway et al., 2008; Finlay et al., 2013), particularly in shallow lakes with moderate P levels where high summer N concentrations stimulate algal growth and cause loss of submerged plants (e.g., Jeppesen et al., 2007; Moss et al., 2013). Indeed, a recent assessment of nutrient sources to Llangorse Lake revealed the importance of reducing N inputs if restoration targets are to be met (May et al., 2010).

The role of climate change in exacerbating the symptoms of eutrophication and confounding recovery efforts cannot be ignored. Climate change is predicted to result in higher water temperatures, shorter periods of ice-cover and longer summer stratification (Jeppesen et al., 2010). Models suggest that lakes with long residence times may experience higher P levels in the future under warmer temperatures (Malmaeus et al., 2006) and shallow lakes may be particularly susceptible. Ecological consequences might include earlier appearance of spring blooming phytoplankton and increased proportions of cyanobacteria. In some systems, negative effects may be compensated by greater predation pressure by zooplankton which is known to be positively temperature dependent. However, fish activity may also increase in warmer temperatures thereby reducing zooplankton populations through increased predation (Moss et al., 2003). In addition, changes in mixing may influence the availability of nutrients in the photic zone and higher temperatures may enhance sediment-P release, whilst higher winter precipitation is likely to enhance nutrient loss from cultivated fields (Battarbee et al., 2008).

An examination of the role of climate change in explaining the shifts in the diatom assemblages of the 12 lakes is beyond the scope of this study. However, detailed studies on two of the lakes, Esthwaite Water (Dong et al., 2012) and Loch Leven (Bennion et al., 2012) have attempted to explore the ways in which nutrients and climate interact on decadal and inter-annual timescales to affect the diatom communities. In these two lakes, the diatom response has been limited despite significant decrease in external nutrient loading. Dong et al. (2012) conclude that while nutrients have been important at Esthwaite Water during the entire 60-year investigation period, air temperature has become a controlling factor in recent decades during a period when nutrient availability was relatively high. Bennion et al. (2012) showed that at an inter-annual scale the diatom data for Loch Leven exhibit high variability, yet there are several changes in species composition in the recent fossil record that may be attributed to climatic controls. In both of these studies the presence of Aulacoseira granulata and Aulacoseira granulata var. angustissma seems to coincide with warmer temperatures. While the diatoms in Lake Mjøsa have experienced shifts toward the pre-enrichment community, not all the pre-eutrophication taxa have reappeared, and analysis of instrumental records lead Hobæk et al. (2012) to conclude that this is either because nitrate concentrations remain high or because water temperature has increased. Such investigations contribute to a better understanding of the effects of multiple environmental drivers on aquatic ecosystems but equally illustrate the complexity of ecosystem response to synchronous changes in nutrients and climate, and the difficulty of disentangling the effects of these interacting pressures (Battarbee et al., 2012). Models that predict likely outcomes of climate change on nutrient regimes will play a vital role in improving our understanding of future lake response and in guiding management decisions (e.g., Whitehead et al., 2006; Trolle et al., 2015). Whilst sediment records cannot be used in a predictive capacity, they provide an opportunity to validate hindcasts derived from dynamic models (Anderson et al., 2006). They should, therefore, play an increasingly important role in assessing uncertainty associated with future predictions.

Conclusions

In terms of the original questions posed we can conclude that the observed changes in the diatom records do reflect both the degradation and the recovery process. The latter has reached a different stage in each of the study lakes and is more clearly seen in several of the deep lakes where the diatom assemblages have started to revert back toward those seen prior to enrichment. In shallow lakes factors such as internal loading and top down control may influence the recovery process and in this study, whilst the assemblages of several shallow lakes were replaced by ones associated with lower productivity following remediation, they did not track back along the enrichment pathway. It can, therefore, be concluded that the deep stratified lakes tend to follow a more predictable recovery pathway than the shallow lakes. Nevertheless, the recovery process has a long way to go in all cases as the present assemblages remain very different from those seen in the pre-enrichment samples and there are indications that for some lakes reference communities may not return. Dissimilarity and principal curve scores are shown to be useful measures for quantifying the deviation from reference condition and provide a useful method with which to evaluate recovery success in the context of global change.

The study highlights the important role that paleolimnological approaches can play in establishing a benchmark against which managers can evaluate the degree to which their restoration efforts are successful. The decadal to centennial timescale adopted here provides the critical temporal context to inform the difficult decisions that emerge for the management of enriched waterbodies. We recognize that this study is based only on diatom responses and our inferences about biological recovery may therefore be biased, especially as diatoms are arguably one of the most sensitive groups in the system and have short response times relative to other assemblages such as macrophytes and fish. Assessments using multiple assemblages are required to evaluate wider ecosystem responses to environmental stressors, hence multi-proxy palaeoecological techniques have an important role to play in future studies of degradation and recovery pathways. The heterogeneous lake dataset used in this study has enabled the identification of broad-scale ecological patterns in recovering temperate lakes. The next step is to refine our understanding of these patterns for different lake types by analysing groups of lakes with more similar attributes.

Author Contributions

HB led the writing of the paper with contributions from GS and BG. HB and BG counted the diatoms in most of the cores and discussed interpretations. GS carried out the numerical analyses. All authors reviewed and edited the manuscript.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank all those colleagues who contributed toward the original projects during which the cores were collected and analyzed. Particular thanks to Peter Appleby and Handong Yang for the radiometric dating, Carl Sayer, Gina Clarke and Samanta Skulmowska for diatom analysis, and to Cath D'Alton of the Cartography Unit, UCL for drawing the figures. We also acknowledge the helpful comments of the reviewers. This work was funded by the EU WISER project (contract no. 226273) and the EU REFRESH project (contract no. 244121). GS was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) Discovery Grant Program (RGPIN 2014–04032).

Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fevo.2015.00094

References

Anderson, N. J. (1997). Historical changes in epilimnetic phosphorus concentrations in six rural lakes in Northern Ireland. Freshwat. Biol. 5, 427–440. doi: 10.1046/j.1365-2427.1997.00249.x