Abstract

Despite the fact that the increased use of elements linked to the Anthropocene is frequently assumed to lead to an increase in the concentrations of the elements in surface waters, temporal trends of trace element (TE) concentrations have rarely been checked. A temporally extended, traceable dataset of TE concentrations in the waters of Lake Geneva, Switzerland (1996–2015) has been used here to explore methodological and data treatment issues that arise when attempting to rigorously determine temporal trends in freshwater TE concentrations. The trace elements studied (Cd, Co, Gd, Mo, Pb, Sb, Sr) have been chosen to cover a wide range of chemical and utilisation conditions. We show that detecting temporal trends from monitoring program data is feasible, even when trends are weak, provided that rigorous data treatment methods are applied. Aspects related to the effect of data quality are discussed in detail. However, ascertaining the statistical significance of any trends calculated remains a difficult issue. With the exception of Co and Sr, that show no significant changes, and Pb, that shows a general decrease, concentrations in lake waters of the trace elements considered have increased significantly, particularly between 2006 and 2015.

INTRODUCTION

Anthropogenic activity has significantly perturbed the Earth's natural biogeochemical cycles, resulting in a significant translocation of many chemical elements from the Earth's crust to human-built objects and to emissions to the environmental compartments. In some environmental media, changes in trace element (TE) concentration have been related to their use. This is the case for sediments, mosses, ombrotrophic peatbogs, or ice cores. It is often assumed that increased human use of chemical elements should also lead to increasing concentrations in surface continental waters (far from direct pollution sources), but this is not clear-cut. Temporal trends in concentrations of major ions have been studied, mostly in the context of acid rain recovery (Garmo et al. 2014; Vuorenmaa et al. 2017) or of organic carbon evolution (Filella & Rodríguez-Murillo 2014), but studies on temporal trends of TE in surface waters are comparatively very scarce (Huser et al. 2011; Todd et al. 2012). The main reason is probably that there are few good, analytically traceable, long-term data series of TE concentrations – like those required for long-term trends to emerge over temporal variations due to occasional changes, seasonality and analytical errors. In particular, TEs have suffered from a lack of appropriate analytical methods until relatively recently (e.g., quadrupole-based inductively coupled plasma–mass spectrometry (ICP-MS) was introduced commercially in 1983) and has benefited from improved trace metal clean procedures introduced in the early 1990s.

Although often just visual (Hatje et al. 2016) or linear correlation methods are applied, rigorous work requires the use of non-parametric methods and the consideration of the possible existence of serial dependence and long-term persistence (LTP). The purpose of the present study is to explore methodological and data treatment issues when attempting to rigorously determine temporal trends in freshwater TE concentrations. For this purpose, we have used a temporally extended, traceable dataset of waters from Lake Geneva, Switzerland. TE concentrations have been collected since 1994 in this lake as part of a water-quality monitoring program set up by the Service de l'Ecologie de l'Eau (SECOE). Lake Geneva is a particularly interesting system because it is a high volume water body with an elevated water residence time, situated in an anthropised environment but not subjected to large direct pollution sources.

Seven elements were chosen for this study (Cd, Co, Gd, Mo, Pb, Sb, Sr). They cover different types of chemical behaviour (Figure SI1) and they also differ with respect to the type and intensity of their use (Figure SI2). In particular, Gd and Pb are interesting reference points. The main use of Pb (i.e., leaded gasoline) was completely phased out in all countries over a relatively short period of time that is well documented (Filella 2017). A noticeable, widespread increase of Gd concentration in surface waters in the past decade has been world-wide observed, likely due to the element's use in magnetic resonance imaging contrast agents (Kulaksiz & Bau 2011; Hatje et al. 2016).

MATERIALS AND METHODS

System, sampling and measurements

Lake Geneva lies at an altitude of 372 m and has a surface area of 580 km2 (Figure SI3). It is the biggest body of freshwater in Western Europe and supplies drinking water for more than 800,000 people. Its average water residence time is 11.3 years. Further details can be found in Table SI1.

Monthly measurement of TE concentrations in point G3 (located at 6.2197° E/46.2994° N; water depth: 72 m) started in December 1995. Measurements are taken at ten depths from 0 to 70 m. Water samples are collected with a PVC Kemmerer sampling bottle lowered in the water from a research vessel. Samples for metal analysis are filtered through 0.45 μm pore-size filters (Whatman® Nuclepore), acidified to pH 1 with Suprapur nitric acid and kept at 2 °C(±2) until analysis. Thus, this study is based on so-called ‘dissolved’ concentrations. For some elements, that are present mostly dissolved in freshwaters, such as antimony (Filella et al. 2009), this is not relevant, but for others, such as lead, that have been shown to be present in variable but never negligible amounts in the particulate fraction (Town & Filella 2002), this fact needs to be taken into consideration when interpreting the results.

Concentrations of TE were measured by ICP-MS on a Fisons PQ2+ before 2006 and on a ThermoFisher X7 II from 2006 to 2015, both in their standard configurations at the Service de l'Ecologie de l'Eau (SECOE) of Canton Geneva, Switzerland. The analytical accuracy was followed by analysing the certified reference materials (CRM) SLRS-4 and SLRS-5 (National Research Council Canada, https://www.nrc-cnrc.gc.ca/eng/solutions/advisory/crm/list_product.html). Detection limits (μg L−1) are: Cd, 0.001; Co, 0.01; Gd, 0.0005; Mo, 0.05; Pb, 0.001; Sb, 0.002; Sr, 0.2. The analytical laboratory is accredited for the analysis of the trace elements in freshwaters (accreditation number: STS245; ISO norm 17025).

Data treatment

Data is shown in Figure SI4. Temporal trends have been evaluated for all data (all depths together; this dataset will be named All data), integrated element load (IEL) and TE concentrations at each depth. IEL is calculated by: (i) dividing the water column into layers, each layer containing all the concentration data measured at equally spaced points; (ii) averaging TE concentrations in each layer to obtain the mass of TE in a column of 1 dm2 of section area and the height of the corresponding layer; and (iii) adding all TE masses from surface to bottom. Integrated TE is then defined as the above mass divided by the lake depth in dm; this represents the weighted TE concentration at the sampling point.

Non-parametric methods are required to analyse temporal trends of TE concentrations in water because data are usually non-normally distributed (Reimann & Filzmoser 2000). Thus, temporal trends have been studied by applying Mann–Kendall (MK) or seasonal Mann–Kendall (SMK) methods (Hirsch & Slack 1984; Helsel & Hirsch 2002), depending on the presence or not of seasonality in the data. Seasonality is expected in lakes, particularly for biogenic elements (Sr, Mo, Co). Here, it has been assessed using the Kruskall–Wallis test of equal medians for the data classified by depth and month. The test was applied to data at 0 and 70 m. Only Sr concentrations at 0 m data showed seasonality (p = 1.5 × 10−6 < 0.05) (p = 0.074 at 70 m depth). Thus, MK was applied for all TE except Sr. IEL data were also assessed for seasonality, which again was found only for Sr (p = 0.019).

The statistical significance of the trends is obtained in MK methods by applying the Z-statistic test of the sum of signs of the differences between every pair of values; Z shows a normal distribution and therefore the statistical significance of the temporal trend can be evaluated by the p value. However, blind use of p values is being increasingly challenged, particularly in the case of tiny effects (Siontis & Ioannidis 2011). We therefore show all trends, whether significant or not, always accompanied by their corresponding Z value, which is directly related to p.

The average change over time of TE concentrations was calculated using Sen's slopes (Sen 1968). A Sen's slope is the median of the slopes between all possible pairs of temporal data (N(N−1)/2), N being the number of pieces of data in the series.

Measured values are very close to detection limits for some elements (Cd, Gd and Pb) and their data series include some censored data (Table SI2). Data have been considered censored when equal to or less than the minimum positive value in each data series. When possible, Sen's slopes have been calculated using the non-parametric method of Akritas et al. (1995) for data at each depth to take the effect of censoring into account.

Median TE concentration profile values have been calculated for the period 2006–2015. Significant differences between Sen's slopes at different depths were assessed by the bootstrap technique following Wilcox (2012). In brief, 10,000 samples (i.e., concentration data series) are obtained by resampling – random picking of elements with replacement – from each time series at each depth, giving 1,000 bootstrap Sen's slopes (SSb). Sen's slope at another depth j (SSj) is then compared with Sen's slope of a depth i (SSi) by determining whether SSj is within 95% of SSb.i. If not, the null hypothesis H0 = ‘equal slopes’ can be rejected at the 1−α confidence level (i.e., SSj and SSi are different). The same technique was applied to estimate the significance of differences between median concentrations at the ten depths.

Hamed and Rao's method (Hamed & Rao 1998) was used to correct MK when there was serial correlation and the application of the method was possible.

Programs for MK and SMK tests are from www.mathworks.com/matlabcentral/fx_files/22389/7/sktt.m, www.mathworks.com/matlabcentral/fileexchange/11190-mann-kendalltau-b-with-sens-method-enhanced/content/ktaub.m adapted by us. Kruskall–Wallis is from PAST (Hammer et al. 2001). Bootstrap calculations have been performed with built-in Matlab functions. Hamed and Rao's method was implemented in Matlab using the ‘Mann-Kendall_Modified_Fix.m’ function, developed by Manjurul Hussain Shourov (https://www.researchgate.net/post/Does_anyone_know_the_package_for_Modified_mann_kendall_trend_test_in_R_software) that corrects a bug in http://www.mathworks.com/matlabcentral/fileexchange/25533-mann-kendall-modified-test. Sen's slope calculations and bootstrap for series with censored data have been calculated using R-3.4.1. For slopes, the ‘cenken’ function in the ‘NADA’ package was used (Helsel 2012).

RESULTS

TE concentrations for Cd, Co, Mo, Pb, Sb, Sr (1995–2015) and Gd (2006–2015) at all depths in Lake Geneva (G3) (All data) are shown in Figure 1. The dramatic improvement in the precision of the measurements due to the change of the ICP-MS instrument in 2006 is readily visible. This effect is noticeable for all TEs considered and is particularly strong for Cd and Sb.

Figure 1

Time series of monthly Cd, Co, Gd, Mo, Pb, Sb and Sr concentrations for 1995–2005 and 2006–2015 in Lake Geneva, Switzerland.

Figure 1

Time series of monthly Cd, Co, Gd, Mo, Pb, Sb and Sr concentrations for 1995–2005 and 2006–2015 in Lake Geneva, Switzerland.

Before any data were treated, the autocorrelation structure of temporal series of IEL data was evaluated. Autocorrelation functions (Figure SI5) show that, except for Cd, significant serial correlation exists for all elements and, consequently, the Hamed and Rao method to evaluate the significance of the temporal trends needs to be applied.

No visually apparent change of tendency is seen along the whole studied period in All data concentrations of three elements (Mo, Sb and Sr) (Figure 2). The change of ICP-MS apparatus in 2006 clearly introduced an offset in the data values for Mo and Sb. Data for these elements were treated separately before and after this year to check for possible spurious effects of the instrument change on the trends observed (Table 1). Although no shift is visible for Sr, data for this element were also treated separately. No significant trends exist for Sr. Increases were observed for All data concentrations and IEL for Mo and Sb in all periods, but they became non-significant for IEL before 2006 when Hamed and Rao's corrected Z values are considered. Gd concentration increases significantly (2006–2015), both considering All data and IEL.

Figure 2

Median concentration profiles and Sen's slopes for Cd, Co, Gd, Mo, Pb, Sb and Sr in Lake Geneva (2006–2015). Error bars in profiles correspond to the interquantile range.

Figure 2

Median concentration profiles and Sen's slopes for Cd, Co, Gd, Mo, Pb, Sb and Sr in Lake Geneva (2006–2015). Error bars in profiles correspond to the interquantile range.

Table 1

Lake Geneva (G3). Results of All data concentration and IEL (integrated element load) trends before and after 2006 (Gd only after 2006). Values in bold are significant (p < 0.05). In IELs, slopes in italics are significant only when Hamed and Rao correction is not applied

  All data Za Slopea/ng L1 y1b IEL Za Hamed Rao Za Slope/μg y1c Yearly increasee/% yr1 
  
Cd         
Before 2006 1,080 0.96 108 −0.16 −0.16 −0.017 (− 0.35–0.27) −0.47 
After 2006 1,150 18.4 0.38 (0.34–0.42) 113 6.70 4.58 0.29 (0.21–0.36) 15 
Co         
Before 2006 1,053 0.11 105 0.094 0.086 0.088 (− 1.4–1.6) 0.14 
After 2006 1,150 7.10 3.1 (2.2–3.9) 114 2.24 1.21 2.2 (0.254.1) 2.6 
Gd         
After 2006 1,093 19.5 0.15 (0.13–0.17) 109 7.26 4.58 0.14 (0.11–0.16) 20 
Mo         
Before 2006 1,080 12.2 30 (25–34) 108 4.24 1.80 21 (11–31) 1.9 
After 2006 1,142 22.1 17 (16–18) 113 7.52 4.24 11 (9.2–13) 1.2 
Pb         
Before 2001 600 −3.67 −1.7 (− 1.7–0) 60 2.23 1.62 4.7 (− 9.9–0.41) 22 
After 2001 1,608 −2.64 −0.0068 (− 0.030–0) 157 0.39 −0.48 −0.11 (− 0.48–0.18) −1.8 
Sb         
Before 2006 1,079 5.83 1.8 (1.2–2.3) 107 1.96 1.21 0.94 (0–1.8) 1.0 
After 2006 1,150 29.6 3.1 (2.9–3.2) 114 10.1 4.60 2.2 (2.0–2.5) 3.7 
Sr         
Before 2006 1,080 1.40 7.5 (7.1–7.7) 108 1.55 d 5.3 (5.255.4) 1.5 
After 2006 1,143 1.38 0.95 (0.93–1.0) 114 1.22 d 0.65 (0.610.67) 0.21 
  All data Za Slopea/ng L1 y1b IEL Za Hamed Rao Za Slope/μg y1c Yearly increasee/% yr1 
  
Cd         
Before 2006 1,080 0.96 108 −0.16 −0.16 −0.017 (− 0.35–0.27) −0.47 
After 2006 1,150 18.4 0.38 (0.34–0.42) 113 6.70 4.58 0.29 (0.21–0.36) 15 
Co         
Before 2006 1,053 0.11 105 0.094 0.086 0.088 (− 1.4–1.6) 0.14 
After 2006 1,150 7.10 3.1 (2.2–3.9) 114 2.24 1.21 2.2 (0.254.1) 2.6 
Gd         
After 2006 1,093 19.5 0.15 (0.13–0.17) 109 7.26 4.58 0.14 (0.11–0.16) 20 
Mo         
Before 2006 1,080 12.2 30 (25–34) 108 4.24 1.80 21 (11–31) 1.9 
After 2006 1,142 22.1 17 (16–18) 113 7.52 4.24 11 (9.2–13) 1.2 
Pb         
Before 2001 600 −3.67 −1.7 (− 1.7–0) 60 2.23 1.62 4.7 (− 9.9–0.41) 22 
After 2001 1,608 −2.64 −0.0068 (− 0.030–0) 157 0.39 −0.48 −0.11 (− 0.48–0.18) −1.8 
Sb         
Before 2006 1,079 5.83 1.8 (1.2–2.3) 107 1.96 1.21 0.94 (0–1.8) 1.0 
After 2006 1,150 29.6 3.1 (2.9–3.2) 114 10.1 4.60 2.2 (2.0–2.5) 3.7 
Sr         
Before 2006 1,080 1.40 7.5 (7.1–7.7) 108 1.55 d 5.3 (5.255.4) 1.5 
After 2006 1,143 1.38 0.95 (0.93–1.0) 114 1.22 d 0.65 (0.610.67) 0.21 

aAll calculated by MK, except Sr (SMK). 95% confidence intervals are for all Sen's slopes' values in MK calculations and for Sen's slopes' seasonal median values in SMK calculations; thus they are not directly comparable.

bExcept for Sr, in μg L−1 y−1.

cExcept for Sr, in mg y−1.

dHamed and Rao correction not needed because SMK includes serial correlation correction.

eYearly increase = 100 × Sen's slope/median IEL (2006).

Changes in direction of the trend during the sampling period are apparently visible for Cd, Co and Pb (Figure 1). Their concentrations appear to decrease in the first half of the period considered and are stationary or increase in the second half. Since the MK method requires a monotonic trend (otherwise the trend could become dependent on the initial and final dates of the time series (Stevenson et al. 2012)), the two data periods were treated separately. Turning years were 2001 for Pb and 2006 for Cd and Co. All data concentrations and IEL of Cd and Co do not change significantly before the turning year and increase significantly afterwards. However, the IEL increasing trend for Co becomes non-significant when autocorrelation is considered. All data concentrations and IEL for Pb decrease in all periods, although the change in IEL is not significant once the Hamed and Rao correction is applied.

Net and percentage changes in IELs can be found in Table 1 for the period 2006–2015.

Lake data at different depths were evaluated only after 2006 in order to avoid noisy data. Median concentration profiles of TEs in the lake (2006–2015) are shown in Figure 2. Results of the bootstrap data treatment are shown in Table SI3. Co, Gd and Mo show no statistically significant changes of their median concentrations with depth while Pb concentrations both at the surface and at the bottom are significantly higher than at other depths. Sb also shows differences (higher concentration) at 0 m compared with 20 and 70 m. Sr median concentration values are higher below 20 m. Cd concentrations are lower at 2.5 m but significantly only when compared with a few other depths (20, 50 and 70 m).

Sen's slopes at all depths (2006–2015) are also shown in Figure 2. Numerical values are given in Table SI4. All are positive, except for Pb at the surface, and statistically significant at all depths in the case of Cd (except at one depth), Gd, Mo and Sb. For Sr, they are only significant at certain depths. Differences in Sen's slope values between depths, evaluated by bootstrap, are shown in Table SI5. Co and Sr show no variation at different depths. Mo, Pb and Sb have lower values at the surface compared to all depths in the case of Pb and to some deeper levels in the case of Mo (7.5, 15, 20 and 30 m) and Sb (15, 20, 50 and 70 m). Some elements show different Sen's slopes at the bottom, smaller in the case of Mo (compared to slopes at 7.5 and 20 m) and bigger for Sb (compared to slopes at 0, 5 and 30 m).

DISCUSSION

Table 2 summarises the results obtained for the various elements, periods and parameters considered. It reflects the complexity of the data treatment, particularly the difficulty of assigning statistical significance to weak trends. The work carried out has made it possible to establish the temporal trends for TEs in Lake Geneva, but also to gain a better understanding of the methodological caveats that apply to temporal trend evaluation of TEs in freshwaters. These are discussed below.

Table 2

Lake Geneva (G3). Temporal trend resume. Results of All data concentration and IEL (integrated element load) trends before and after 2006 (excep Pb, before and after 2001, and Gd only after 2006) and Sen's slopes after 2006. Values in bold are significant (p < 0.05). Slopes in italics are significant only when Hamed–Rao correction is not applied. Remember that this correction could not be applied to All data

Element All data
 
IEL
 
Sen slopes at different depths 
First period Second period First period Second period >2006 only 
Cd Increase Decrease Increase Increasea 
Co Increase Increase Increase Increaseb 
Gd – Increase – Increase Increase 
Mo Increase Increase Increase Increase Increase 
Pb Decrease Decrease Decrease Decrease 0 m: decrease; rest: increase 
Sb Increase Increase Increase Increase Increase 
Sr Increase Increase Increase Increase Increase 
Element All data
 
IEL
 
Sen slopes at different depths 
First period Second period First period Second period >2006 only 
Cd Increase Decrease Increase Increasea 
Co Increase Increase Increase Increaseb 
Gd – Increase – Increase Increase 
Mo Increase Increase Increase Increase Increase 
Pb Decrease Decrease Decrease Decrease 0 m: decrease; rest: increase 
Sb Increase Increase Increase Increase Increase 
Sr Increase Increase Increase Increase Increase 

aNot statistically significant at 20 m depth.

bStatistically significant at 10, 20, 30 and 50 m depth.

Methodological issues

Need for data traceability

Data in trend studies have usually been obtained from monitoring programs rather than from research-oriented projects, which implies different objectives and data quality standards. Nevertheless, our study shows that temporal trends can be followed using monitoring data, provided that it is possible to have direct access to the raw data and protocols used (i.e., changes of instrument, calibration, use of CRMs, data pre-treatment). However, detailed information on quality assurance (QA) all through the years is not available even in a ‘favourable’ case like ours. Unfortunately, this is a problem inherent to any data trend study based on monitoring data.

Data quality and length of the data series

How long does a data series need to be for a temporal trend to be detected? This depends on two factors: the strength of the change and the quality of the data. It can be roughly evaluated by calculating the time of emergence (ToE) parameter, which gives the point in time when a signal finally emerges from the background noise of natural variability (Keller et al. 2014). Values for ToE of the elements studied (Table 3), calculated before and after 2006, clearly show that where data are noisier (i.e., before 2006), longer time periods are needed for the trend to be detectable. In the case of Pb and Sr, ToEs are higher in the second period because trends are stronger in the first period and compensate for the greater noise. As expected, elements with significant changes have shorter ToEs, closer to ten years.

Table 3

Times of emergence (ToE) of IELs

Element ToE/ya
 
< 2006 > 2006 
Cd 1,373 13.6 
Co 822 42.5 
Gd  12.9 
Mo 22.8 20.5 
Pb 63.5 744 
Sb 55.5 7.67 
Sr 30.0 75.5 
Element ToE/ya
 
< 2006 > 2006 
Cd 1,373 13.6 
Co 822 42.5 
Gd  12.9 
Mo 22.8 20.5 
Pb 63.5 744 
Sb 55.5 7.67 
Sr 30.0 75.5 

aToE is defined as ToE = N/S where N is a measure for variability (in our case, percentile range 0.025–97.5, in μg) and S the trend (Sen's slope in μg y−1) (Keller et al. 2014).

Effect of censoring

Pb, Cd and Gd data have a sizeable fraction of censored data (Table SI2). Simply substituting a fabricated value (zero, detection limits, or another value) for censored data can lead to wrong results in statistical data (mean, standard deviation), hypothesis tests and regression models (Helsel 2012). Median concentrations are not affected by left-censored data (i.e., those less than a limit) if censoring is less than 50% of the data. Although trends are most effectively detected in uncensored data as compared to censored data, the existence of censored data will not give any spurious significant trend (i.e., its effect is eliminating information). This means that statistical significance in MK is not affected by censoring. However, it has an effect on the value of Sen's slopes (Table SI2). Those slopes have thus been calculated using the non-parametric method of Akritas et al. (1995) for data at each depth, but calculations could not be performed for All data because the method does not work when tied values are present (in our case, there are ten concentration values, i.e., ties at each time).

Other issues related to data quality

Another effect of data quality is related to the precision of the measurements that, when close to the detection limit, results in very few significant figures in the readings. Here lies the reason behind the zero Z values for Cd and Co data before 2006 (remember that MK is a ranking method) and the ‘ugly’ aspect of Gd in Figure 1. A final aspect that has, to our knowledge, never been explored in concentration trend studies, is the effect on the measured concentration values of potential drifts of the instruments over the years. We have tried to follow potential drifts by applying the same data treatment to the CRM concentrations measured (Table SI6), with no conclusive results (see Discussion in the Supplementary material).

Serial correlation and persistence

Finding statistically significant, deterministic temporal trends in time series is complicated by two features of time series: serial correlation and LTP. Hydrological time series, such as those we are dealing with, very often show positive serial correlation. Serially correlated data without any deterministic trend tend to stick above or below the average value of the series more frequently than would be the case for a purely random process. This increases the likelihood of identifying a significant trend where there is none. Hamed and Rao's method (Hamed & Rao 1998) was used to correct MK when there was serial correlation, but it could only be applied to IEL data and not to data from all depths lumped together because the method does not allow for ties (more than one value at each time). As expected (Cox & Stuart 1955), the application of the corresponding correction unequivocally leads to a decrease in the significance of the observed trends when it could be applied (i.e., IEL).

A ‘persistent’ time series shows long periods of values above or below the mean, even if the mean does not change over the long term (Beran 1994). Unfortunately, deciding whether LTP is present and accounting for it requires data series longer than the ones we are studying here (O'Connell et al. 2016). We therefore did not correct our results for LTP. A deeper discussion on these problems can be found in the SI file (Section SI1).

Looking for causes

Trend studies commonly attempt to ascribe observed effects to well-defined, physical causes. This is usually done by establishing correlations, and backing them up by informally including other lines of evidence from the literature in the discussion section (Downes et al. 2002). Apart from the well-known fact that correlation is not causation, in practice, this approach is not straightforward. Only in rare cases is the trend strong enough and due to a few, well-defined, drivers (e.g., recovery from acid rain); generally, observed effects stem from the interaction of a variety of interrelated causes that are not always acting in the same direction, often leading to ‘weak’ causation effects. This type of exercise or the application of more sophisticated approaches (Granger 1980; Kleinberg 2012) is entirely outside of the scope of this paper. It is interesting however to discuss what effects can be reasonably expected and whether the results obtained fit with existing knowledge.

When an input will have a measurable effect

Temporal changes in lake TE concentrations will depend on temporal changes in lake TE inputs (i.e., direct input from waste water treatment plants, direct atmospheric deposition on the lake, input from the lake catchment) and/or temporal changes of internal physical and chemical processes taking place in the water column and sediments. Obviously, not every change in an element input will result in a statistically significant increase in concentrations because, on the one hand, increasing input load might be compensated for by internal lake processes while, on the other, the input itself must be (i) big enough compared to other sources, (ii) big enough compared to lake concentrations and (iii) present in a ‘convenient’ form. Examples of these conditions are found in our results.

Big enough compared to other sources: There are 171 WWTPs in the Lake Geneva watershed, which is home to a 3,009,830 population equivalent (Table SI2). A recent quantitative assessment of TE flows from WWTP in Switzerland (Vriens et al. 2017) showed that the contribution of WWTP effluents to river flow is very low for most of the elements we are concerned with here: Cd (3%), Co (8%), Mo (5%), Pb (1%) and Sr (2%). The WWTP contribution is a bit higher for Sb (16%) and huge for Gd (83%). In the case of Gd, it is widely accepted that this is due to its medical use (Kuroda et al. 2016) and, for Sb, it is probably the result of its widespread use and dissipation in the final products (Turner & Filella 2017). Therefore, the only case where the observed temporal trend can be unequivocally attributed to the input from WWTP is Gd.

Big enough compared to lake concentrations: Direct deposition of TE on the lake will have a measurable effect on surface concentrations only when deposition flux and accumulation thereof are at least an order of magnitude higher than element water concentrations. For instance, if we consider the elements studied here, Sr concentrations in the lake are more than two to five orders of magnitude higher than for the other elements studied, while atmospheric deposition is negligible. It is therefore only to be expected that no effect is observed. Atmospheric deposition would need to be very high indeed for an effect to be seen for this element. In fact, the only element that shows statistically significant higher concentrations at the surface compared to other depths (2006–2015) is Pb, directly pointing to the importance of atmospheric deposition for this element. Atmospheric deposition cannot be excluded for the other elements, but it is not equally translated into higher surface concentrations, except for Sb, albeit less clearly than for Pb. Atmospheric input is probably not ‘seen’ for Mo because of its higher concentration in lake water compared to Pb and Sb, although Mo concentrations in ice cores in Mount Rosa (Swiss–Italian Alps) – a good measure of metal deposition – were similar to Sb ones (at least before 1995) (Barbante et al. 2004).

Present in a ‘convenient’ form: According to these relatively old data in Mount Rosa, Co deposition was about five times Sb deposition and, in spite of both elements having similar concentrations, no change is observed in Co. This could be explained by the higher variability of Co data (Figure 1) or by Co being present mostly in not readily soluble dust, mainly of Saharian origin (Heimbürger et al. 2010). When TE deposition is in aerosol or particulate form, any effect on lake dissolved concentrations will require that the element is readily, or at least partially, dissolved in lake water. This might be limiting not only in the case of Co, but also for Pb in re-suspension of fine lead-contaminated soil dust (Laidlaw et al. 2012) or of the Sb released from brake pads (Iijima et al. 2008).

Further information from the treatment of lake layers

Temporal changes of TE concentrations in the surface and upper lake layers, as expressed in Sen's slopes, provide additional information on changes in atmospheric deposition. The most extreme effect on surface Sen's slopes is observed in Pb, where the slope changes to a negative, although non-significant, value. Sb and Mo also show (significant) smaller surface slopes. The simplest explanation of smaller Sen's slopes in and near the surface is the decrease of atmospheric deposition on the lake during the period considered. Detailed changes of TE atmospheric deposition on Lake Geneva during 1995–2015 are not known but existing data for Switzerland point to a general decrease. For instance, drastic decreases of Pb (86%) and Cd (66%) deposition and slight decreases for Co have been measured in Swiss mosses (1990–2010) (Harmens et al. 2013). A 23% average decrease was observed for Sb in mosses from seven European countries from 2005 to 2010 (Harmens et al. 2013).

Discrepant results and some qualitative speculation

Even if current decreasing atmospheric TE deposition following emission control is widely accepted (Harmens et al. 2010), except for Pb, our results point to an increase, not a decrease in TE lake concentrations. Thus, a further source seems to be at play. A plausible candidate is a change in the release from the lake catchment.

Not all elements are equally released from soils to waters. An indication of their soil mobility is given by the ratio between their concentrations in soils and in waters. These ratios, estimated from median concentration values in European soils (FOREGS) (Salminen et al. 2005) (soil concentrations in mg g−1 and water concentrations in mg L−1), give: Gd (385) > Pb (243) > Co (49) > Cd (15) > Sb (8.6) > Mo (2.8) > Sr (0.81). Very low mobility is thus to be expected for ‘natural’ Gd and Pb as compared to the other elements. Release of TEs deposited from atmospheric deposition or directly added to the soils (e.g., the case of Cd introduced by the use of fertilisers (Roberts 2014)) might differ depending on the form of the deposited element and its interaction with soils but similar behaviour to ‘natural’ TEs is likely. For instance, it is well known that Pb in soils coming from its use in leaded gasoline is extremely insoluble and immobile, with an expected ‘half-life’ of 700 years (Semlali et al. 2004). In any case, changes in TE release from soils, originating either from factors that affect rock weathering (effect of climate change) (Taylor et al. 2012) or from release of TE deposited from atmospheric deposition, are very slow and a delay of years or decades can be expected (Lawlor & Tipping 2003; Huser et al. 2011). Therefore, the origin of the increasing TE concentrations in water (Cd, Mo and Sb) lies very probably in the release of TE previously deposited in soils. A similar effect has been detected in streams and rivers of southern Sweden for Pb (Huser et al. 2011).

CONCLUSIONS

Rigorous evaluation of temporal trends of TE in surface waters is not trivial because good, long, traceable concentration data series are required to evaluate whether trends exist. Unfortunately, these data series are scarce, are relatively short because the necessary analytical techniques are quite new, and the quality of the data has until recently often been poor. Moreover, these data have usually been obtained under monitoring programs, not from research-oriented projects, which implies different data quality standards. Nevertheless, studying temporal trends is feasible, provided direct access to raw data is possible and rigorous data treatment methods are applied.

Gd constitutes a good example of how direct input into environmental compartments may have a direct impact, which is easy to detect and explain. However, this type of situation is more the exception than the rule. Neither the often assumed general increase of TE concentrations in environmental compartments linked to the current massive use of the elements, nor the expected decrease, as a result of general, diminishing atmospheric deposition of TE in Europe, are confirmed by this study. This is a nice illustration of the complexity of the response of natural systems to disruption and a warning against over-quick generalisations in environmental sciences.

REFERENCES

REFERENCES
Akritas
,
M. G.
,
Murphy
,
S. A.
&
LaValley
,
M. P.
1995
The Theil-Sen estimator with doubly censored data and applications to astronomy
.
Journal of the American Statistical Association
90
,
170
177
.
Barbante
,
C.
,
Schwikowski
,
M.
,
Döring
,
T.
,
Gägeler
,
H. W.
,
Schotterer
,
U.
,
Tobler
,
L.
,
Van De Velde
,
K.
,
Ferrari
,
C.
,
Cozzi
,
G.
,
Turetta
,
A.
,
Rosman
,
K.
,
Bolshov
,
M.
,
Capodaglio
,
G.
,
Cescon
,
P.
&
Boutron
,
C.
2004
Historical record of European emissions of heavy metals to the atmosphere since the 1650s from Alpine Snow/Ice Cores Drilled near Monte Rosa
.
Environmental Science & Technology
38
,
4085
4090
.
Beran
,
J.
1994
Statistics for Long-Memory Processes
.
CRC Press
,
Boca Raton, FL
,
USA
, pp.
315
.
Downes
,
B. J. L.
,
Barmuta
,
A.
,
Fairweather
,
P. G.
,
Faith
,
D. P.
,
Keough
,
M. J.
,
Lake
,
P. S.
,
Mapstone
,
B. D.
&
Quinn
,
G. P.
2002
Monitoring Ecological Impacts: Concepts and Practice in Flowing Waters
.
The University of Chicago Press
,
Cambridge
,
UK
.
Filella
,
M.
2017
Environmental impact of alkyl lead(IV) derivatives: perspective after their phase-out
.
Metal Ions in Life Sciences
17
,
471
490
.
Filella
,
M.
,
Williams
,
P. A.
&
Belzile
,
N.
2009
Antimony in the environment: knowns and unknowns
.
Environmental Chemistry
6
,
95
105
.
Garmo
,
Ø. A.
,
Skjelkvåle
,
B. L.
,
de Wit
,
H. A.
,
Colombo
,
L.
,
Curtis
,
C.
,
Fölster
,
J.
,
Hoffmann
,
A.
,
Hruška
,
J.
,
Høgåsen
,
T.
,
Jeffries
,
D. S.
,
Keller
,
W. B.
,
Krám
,
P.
,
Majer
,
V.
,
Monteith
,
D. T.
,
Paterson
,
A. M.
,
Rogora
,
M.
,
Rzychon
,
D.
,
Steingruber
,
S.
,
Stoddard
,
J. L.
,
Vuorenmaa
,
J.
&
Worsztynowicz
,
A.
2014
Trends in surface water chemistry in acidified areas in Europe and North America from 1990 to 2008
.
Water, Air, & Soil Pollution
225
,
1880
.
DOI 10.1007/s11270-014-1880-6
.
Granger
,
C. W. J.
1980
Testing for causality: a personal viewpoint
.
Journal of Economic Dynamics and Control
2
,
329
352
.
Hamed
,
K. H.
&
Rao
,
A. R.
1998
A modified Mann–Kendall trend test for autocorrelated data
.
Journal of Hydrology
204
,
182
196
.
Hammer
,
Ø.
,
Harper
,
D. A. T.
&
Ryan
,
P. D.
2001
PAST: Paleontological Statistics Software Package for Education and Data Analysis
.
Palaeontologia Electronica
, pp.
4
.
Harmens
,
H.
,
Norris
,
D. A.
,
Steinnes
,
E.
,
Kubin
,
E.
,
Piispanen
,
J.
,
Alber
,
R.
,
Aleksiayenak
,
Y.
,
Blum
,
O.
,
Coşkun
,
M.
,
Dam
,
M.
,
De Temmerman
,
L.
,
Fernández
,
J. A.
,
Frolova
,
M.
,
Frontasyeva
,
M.
,
González-Miqueo
,
L.
,
Grodzińska
,
K.
,
Jeran
,
Z.
,
Korzekwa
,
S.
,
Krmar
,
M.
,
Kvietkus
,
K.
,
Leblond
,
S.
,
Liiv
,
S.
,
Magnússon
,
S. H.
,
Mankovská
,
B.
,
Pesch
,
R.
,
Rühling
,
A.
,
Santamaria
,
J. M.
,
Schröder
,
W.
,
Spiric
,
Z.
,
Suchara
,
I.
,
Thöni
,
L.
,
Urumov
,
V.
,
Yurukova
,
L.
&
Zechmeister
,
H. G.
2010
Mosses as biomonitors of atmospheric heavy metal deposition: spatial patterns and temporal trends in Europe
.
Environmental Pollution
158
,
3144
3156
.
Harmens
,
H.
,
Norris
,
D.
&
Mills
,
G.
,
the participants of the Moss Survey
2013
Heavy Metals and Nitrogen in Mosses: Spatial Patterns in 2010/2011 and Long-Term Temporal Trends in Europe
.
ICP Vegetation Programme Coordination Centre, Centre for Ecology and Hydrology
,
Bangor
,
UK
, pp.
63
.
Heimbürger
,
L.-E.
,
Migon
,
C.
,
Dufour
,
A.
,
Chiffoleau
,
J.-F.
&
Cossa
,
D.
2010
Trace metal concentrations in the North-western Mediterranean atmospheric aerosol between 1986 and 2008: seasonal patterns and decadal trends
.
Science of the Total Environment
408
,
2629
2638
.
Helsel
,
D. R.
2012
Statistics for Censored Environmental Data Using Minitab_and R
,
2nd edn
.
Wiley
,
Hoboken, NJ
, pp.
324
.
Helsel
,
D. R.
&
Hirsch
,
R. M.
2002
Statistical Methods in Water Resources; Techniques of Water-Resources Investigations of the United States Geological Survey
.
Book 4, Chapter A3
,
USGS
,
Reston, VA
.
Hirsch
,
R. M.
&
Slack
,
J. R.
1984
A nonparametric trend test for seasonal data with serial dependence
.
Water Resources Research
20
,
727
732
.
Huser
,
B. J.
,
Köhler
,
S. J.
,
Wilander
,
A.
,
Johansson
,
K.
&
Fölster
,
J.
2011
Temporal and spatial trends for trace metals in streams and rivers accross Sweden (1996–2009)
.
Biogeosciences
8
,
1813
1823
.
Iijima
,
A.
,
Sato
,
K.
,
Yano
,
K.
,
Kato
,
M.
,
Kozawa
,
K.
&
Furuta
,
N.
2008
Emission factor for antimony in brake abrasion dusts as one of the major atmospheric antimony sources
.
Environmental Science & Technology
42
,
2937
2942
.
Keller
,
K. M.
,
Joos
,
F.
&
Raible
,
C. C.
2014
Time of emergence of trends in ocean biogeochemistry
.
Biogeosciences
11
,
3647
3659
.
Kleinberg
,
S.
2012
Causality, Probability, and Time
.
Cambridge University Press
,
New York
.
Kuroda
,
K.
,
Itten
,
R.
,
Kovalova
,
L.
,
Ort
,
C.
,
Weissbrodt
,
D. G.
&
McArdell
,
C. S.
2016
Hospital-use pharmaceuticals in Swiss waters modeled at high spatial resolution
.
Environmental Science & Technology
50
,
4742
4751
.
Laidlaw
,
M. A. S.
,
Zahran
,
S.
,
Mielke
,
H. W.
,
Taylor
,
M. P.
&
Filippelli
,
G. M.
2012
Re-suspension of lead contaminated urban soil as a dominant source of atmospheric lead in Birmingham, Chicago, Detroit and Pittsburgh, USA
.
Atmospheric Environment
49
,
302
310
.
O'Connell
,
P. E.
,
Koutsoyiannis
,
D.
,
Lins
,
H. F.
,
Markonis
,
Y.
,
Montanari
,
A.
&
Cohn
,
T.
2016
The scientific legacy of Harold Edwin Hurst (1880–1978)
.
Hydrological Sciences Journal
61
,
1571
1590
.
Salminen
,
R.
,
Batista
,
M. J.
,
Bidovec
,
M.
,
Demetriades
,
A.
,
De Vivo
,
B.
,
De Vos
,
W.
,
Duris
,
M.
,
Gilucis
,
A.
,
Gregorauskiene
,
V.
,
Halamic
,
J.
,
Heitzmann
,
P.
,
Lima
,
A.
,
Jordan
,
G.
,
Klaver
,
G.
,
Klein
,
P.
,
Lis
,
J.
,
Locutura
,
J.
,
Marsina
,
K.
,
Mazreku
,
A.
,
O'Connor
,
P. J.
,
Olsson
,
S. Å.
,
Ottesen
,
R.-T.
,
Petersell
,
V.
,
Plant
,
J. A.
,
Reeder
,
S.
,
Salpeteur
,
I.
,
Sandström
,
H.
,
Siewers
,
U.
,
Steenfelt
,
A.
&
Tarvainen
,
T.
2005
Geochemical Atlas of Europe. Part 1: Background Information, Methodology and Maps
.
Geological Survey of Finland
,
Espoo
, pp.
526
.
Semlali
,
R. M.
,
Dessogne
,
J.-B.
,
Monna
,
F.
,
Bolte
,
J.
,
Azimi
,
S.
,
Navarro
,
N.
,
Denaix
,
L.
,
Loubet
,
M.
,
Chateau
,
C.
&
van Oort
,
F.
2004
Modeling lead input and output in soils using lead isotopic geochemistry
.
Environmental Science & Technology
38
,
1513
1521
.
Sen
,
P. K.
1968
Estimates of regression coefficient based on Kendall's tau
.
Journal of the American Statistical Association
63
,
1379
1389
.
Siontis
,
G. C. M.
&
Ioannidis
,
J. P. A.
2011
Risk factors and interventions with statistically significant tiny effects
.
International Journal of Epidemiology
40
,
1292
1307
.
Stevenson
,
K.
,
Alessa
,
L.
,
Altaweel
,
M.
,
Kliskey
,
A. D.
&
Krieger
,
K. E.
2012
Minding our methods: how choice of time series, reference dates, and statistical approach can influence the representation of temperature change
.
Environmental Science & Technology
46
,
7435
7441
.
Taylor
,
L. L.
,
Banwart
,
S. A.
,
Valdes
,
P. J.
,
Leake
,
J. R.
&
Beerling
,
D. J.
2012
Evaluating the effects of terrestrial ecosystems, climate and carbon dioxide on weathering over geological time: a global-scale process-based approach
.
Philosophical Transactions of the Royal Society B: Biological Sciences
367
,
565
582
.
Todd
,
A. S.
,
Manning
,
A. H.
,
Verplanck
,
P. L.
,
Crouch
,
C.
,
McKnight
,
D. M.
&
Dunham
,
R.
2012
Climate-change-driven deterioration of water quality in a mineralized watershed
.
Environmental Science & Technology
46
,
9324
9332
.
Town
,
R. M.
&
Filella
,
M.
2002
Size fractionation of trace metals in freshwaters: implications for understanding their speciation and fate
.
Re/Views in Environmental Science and Bio/Technology
1
,
277
297
.
Turner
,
A.
&
Filella
,
M.
2017
Field-portable-XRF reveals the ubiquity of antimony in plastic consumer products
.
Science of the Total Environment
584–585
,
982
989
.
Vriens
,
B.
,
Voegelin
,
A.
,
Hug
,
S. J.
,
Kaegi
,
R.
,
Winkel
,
L. H. E.
,
Buser
,
A. M.
&
Berg
,
M.
2017
Quantification of element fluxes in wastewaters: a nationwide survey in Switzerland
.
Environmental Science & Technology
51
,
10943
10953
.
Vuorenmaa
,
J.
,
Augustaitis
,
A.
,
Beudert
,
B.
,
Clarke
,
N.
,
de Wit
,
H. A.
,
Dirnböck
,
T.
,
Frey
,
J.
,
Forsius
,
M.
,
Indrikson
,
I.
,
Kleemola
,
S.
,
Kobler
,
J.
,
Krám
,
P.
,
Lindroos
,
A.-J.
,
Lundin
,
L.
,
Ruoho-Airola
,
T.
,
Ukonmaanaho
,
L.
&
Váňa
,
M.
2017
Long-term sulphate and inorganic nitrogen mass balance budgets in European ICP Integrated Monitoring catchments (1990–2012)
.
Ecological Indicators
76
,
15
29
.
Wilcox
,
R.
2012
Robust regression
. In:
Statistical Modeling and Decision Science: Introduction to Robust Estimation and Hypothesis Testing
,
3rd edn
.
Elsevier
,
Amsterdam
, pp.
471
530
.

Supplementary data