Spatiotemporal variabilities of major, minor and trace elements in the German part of the river Rhine

Abstract

Recently developed multi-element methods for filtered and unfiltered river water enable the analysis of 67 elements in a single analytical run, facilitating the assessment of multi-element fingerprints in monitoring programs. To elucidate the occurrence and possible pathways of emerging contaminants, it is essential to have baseline data for as many elements as possible to detect anomalies at specific environmental events. However, the variability of element concentrations in river water on the temporal and spatial scale is only poorly understood, causing considerable uncertainty for river water monitoring. Therefore, we conducted comprehensive sampling campaigns to assess the spatiotemporal variabilities of element concentrations in the river Rhine (Germany) at different discharge levels. Both, the long-term temporal (one-year data) and the spatial variability data revealed a distinct behavior of two element groups in relation to the discharge: elements showing a dilution effect (e.g.: B, Mg, S, K, Ca, Br, Sr, Mo, and U) or a co-rising effect (e.g.: Al, Si, P, Ti, V, Mn, Fe, Ni, Cu, Ga, As, Rb, Y, Cs, Ba, La, Ce, Pr, Nd, Sm, Gd, Pb, and Th) with higher discharge in the unfiltered river water. However, certain elements such as K, Ba or Gd displayed a variable behavior throughout the dataset, underlining the importance of collecting enough baseline data from various locations and conditions to detect anomalies. The analysis of unfiltered water samples in comparison with the filtered fraction allowed us to detect opposite behavior between the dissolved and particulate fractions for Co, Ni, As, Rb, Cs, and Gd. Cross-profile measurements were conducted to investigate spatial variabilities, revealing that lateral spatial gradients in element concentrations (up to a factor of 4) are more pronounced than depth gradients, likely caused by insufficient mixing of river influents or point-sources. Thus, fixed monitoring stations or single-point sampling for long-term data acquisition might not be able to capture the whole picture regarding element behavior in rivers.

---

Environmental significance

Present challenges in describing the chemical status of rivers are a precise and fast response to rapid changes and the ability to conduct timely forecasts. Therefore, adequate validated multi-analyte methods and sampling procedures must be in place. In previous publications, we presented versatile multi-element methods for (un)filtered river water samples. Here, these methods are applied to provide an in-depth investigation of spatiotemporal variabilities at different locations of the river Rhine at various discharge levels. In summary, low-resolution sampling in time and space disregards the impact of critical point-sources, hides information and thereby induces high uncertainties in monitoring datasets. The multi-element approaches enable the creation of baseline datasets for each element to detect anomalies and to serve as initial training data for forecasting tools.

1 Introduction

Rivers are transboundary vectors for various substances – from inland waters to oceans, across country borders, or from spring- and wastewater to drinking water. The chemical anthropogenic impact on surface waters cannot always be foreseen due to a lack of monitoring data. Moreover, we are currently dealing with contaminations of the past and present and are trying to address future challenges.^1,2 In order to operate an effective alert system for river water quality, long-term baseline data and stable sample processing systems need to be established. It is important to improve our understanding of the transport mechanisms and pathways of substances in rivers and to have monitoring and prediction tools in place, to maintain a healthy environment in rural and industrialized areas. Therefore, it is necessary to be able to distinguish between anthropogenic inputs and near-natural element cycles, as well as to understand the uncertainties of existing monitoring approaches. Separating the different sources contributing to river water chemistry and providing artificial intelligence (AI) training datasets for future data analysis as well as forecasts are among the reasons to produce reliable monitoring data of as many elements as possible at high resolution.^3–6

For inorganic analytes, mainly Cr, Ni, Cu, Zn, As, Cd, Sb, Hg and Pb are included in long-term monitoring databases.^7–9 In addition to these so-called “legacy elements”, rare earth elements (REEs) are also frequently used in multiple application fields and have been found to enter aquatic systems through several pathways, e.g. medical applications, mining, electronic waste or chemical fertilizers.^10–14 However, the number of studies providing REE data from rivers and other aquatic systems is limited.¹⁴ REEs are included in the list of so-called “technology-critical elements” (TCEs) which also encompasses the rarely analyzed elements Ga, Ge, In, Te, Nb, and Ta. They are expected to be found in higher concentrations due their use in modern technology and the resulting increased mining activities.⁷ Furthermore, changing element compositions of rivers can also provide valuable information on near-natural processes within the catchment.¹⁵

In previous studies, we developed multi-element methods for filtered river water samples (<0.45 µm)¹⁶ and whole water samples (unfiltered river water)¹⁷ as tools for a better understanding of the comprehensive element content in rivers (67 different elements). These multi-element methods were developed to be applied in long-term monitoring approaches and for the detection of anomalies over time. Whole water samples enable the assessment of water status including the dissolved, colloidal, and particulate fractions in one sample type. At low suspended particulate matter (SPM) concentrations, when sufficient SPM collection is challenging, the whole water sample is a good alternative, by subtracting the filtered element concentration from the whole water element concentration. Additionally, the whole water sample is the only sample type allowing us to calculate river element budgets including all fractions. Moreover, the multi-element methods were developed to gain knowledge about spatiotemporal variabilities of element concentrations in riverine systems and to improve our understanding of element cycles.

For multiple (most often economical and logistical) reasons, the analysis of spatial or temporal trends along a river is usually undertaken with single-point sampling at a defined depth.^14,18–22 However, it is well-known that SPM and elemental concentrations in rivers can be highly heterogeneous over time and across river cross-sections. For instance, Bouchez et al.²³ conducted several sampling campaigns covering multiple depth profiles in the Amazon River basin. An underestimation by a factor of 2 was found for Al and Cs in the particulate fraction of surface samples compared to depth integrated values.²⁴ To date, sampling campaigns covering several sampling points in the river cross profile for elemental analysis have mostly been conducted for large rivers such as the Amazon River,^23–25 the Ganga basin,²⁶ the Athabasca River,²⁷ or several of the large Arctic Rivers.^28,29 Nevertheless, the availability of long-term data from a single sampling point over a long period of time is also valuable to be able to detect changes over decades. In this case, changing the sampling location is not recommended.⁹ The frequency of sampling is an important question regarding the necessary resolution to capture temporal variabilities.^30–32 Floury et al.³⁰ compared 40 min sampling frequencies with 7 h and daily frequencies using ion chromatography and an automated monitoring station next to the Avenelles River which is part of the Orgeval catchment in France. The results for Ca, as one example, indicate that the lowest Ca concentration was observed during the highest water discharge at rain events by sampling at 40 min intervals. The 7 h and daily sampling frequencies would not have captured this. Therefore, low-frequency sampling could lead to a loss on information.^30–32

In this study, large and small-scale spatiotemporal variabilities in riverine systems were addressed by conducting four different types of sampling campaigns in the German part of the Rhine:

1. One-year temporal variability (daily integrated samples of the first seven days of each month between May 2021 and April 2022).

2. Short-term temporal variability (1-minute resolution for 10 min).

3. Spatial variability (5 vertical samples taken simultaneously at 50 cm intervals across the full river cross-section including lateral and depth variability).

4. Small-scale spatial variability (a 4 × 4 sampling grid covering a 1.6 × 1.6 m window, sampling from one riverbank to the other at a defined depth).

The sampling campaigns were conducted at three different locations along the Rhine, with different hydrological, environmental and anthropogenic conditions and at different discharge stages (low, moderate and high). The Rhine is, at least in comparison to previous literature on spatial or spatiotemporal variability campaigns, a shallow river which might lead to a different effect regarding the river size.^23–27 The main goals of our study include the better understanding (i) of spatiotemporal variabilities in element concentrations in our monitoring approaches and (ii) of similar or differing element behavior in relation to discharge and/or SPM. Samples from each sampling campaign were analyzed for 67 elements, covering major, minor and trace elements, in two different fractions: the filtered (<0.45 µm) and the whole water sample (unfiltered sample). In the following, the term “dissolved” is used for the <0.45 µm fraction next to “filtered”. For the one-year temporal variability study, daily integrated samples were collected over a one-year period between May 2021 and April 2022. To reduce the total number of samples, only the first seven days of each month were investigated. The short-term temporal sampling campaigns were conducted with a sampling frequency of 1-minute intervals over a total time of 10 min. For the spatial variability campaigns, two different measuring frames were developed to enable in situ and simultaneous sampling of several samples by an automatic closing mechanism. One measuring frame was able to hold 5 vertically fixed bottles at 50 cm distances, as well as sensors to determine the water depth. It was used to determine the cross-sectional variability. The second frame was built as a 4 × 4 raster holding 16 bottles at 40 cm distances to cover small-scale spatial variabilities.

Consequently, the different sampling campaigns addressed four different topics. Sampling type no. 1 (one-year data) was used to establish a grouping of element behavior based on samples from regular monitoring networks. Sampling type no. 2 (short-term temporal variability) was conducted to assess the influence of temporal effects on manual consecutive sampling at the same sampling spot and at a defined depth, which is often conducted for grab sampling or generally for single-point sampling. The objective was to examine whether observed spatial variabilities can be to some extent an expression of temporal variabilities. The spatial variability campaigns (no. 3) were used to extend our knowledge regarding multi-point vs. single-point sampling to re-evaluate current monitoring approaches which mostly rely on single-point sampling. The small-scale spatial variability campaign (no. 4) examines the required sampling density to resolve pronounced variabilities in the river cross-profile. To the authors' knowledge, there are no existing studies with similarly dense spatiotemporal sampling campaigns in the river cross-profile including the analysis of 67 major, minor and trace elements in the filtered and unfiltered water fractions. Moreover, the investigation of turbulence-driven small-scale temporal and small-scale spatial variability provided novel comprehensive datasets which are not yet available in the present literature.

2 Materials and methods

2.1 Sampling locations

The Rhine is a major European river with a length of about 1250 km between its source in Switzerland and the North Sea in the Netherlands (of which 870 km flow through Germany). The Rhine catchment area is known to contain densely populated, farming and industrialized sections, leading to a diverse mix of substances originating from different sources.^9,33,34 In this study, all sampling sites were located along the German part of the Rhine (refer to Fig. 1: Koblenz (KAA), Niederlahnstein (LST), Brohl (BRL) and Emmerich (EMM)). They were chosen to cover the different hydrological conditions in the middle and lower Rhine areas as well as the direct influence of the river Moselle after its confluence at the BRL station. Moreover, the site selection was based on available sampling infrastructure from the Federal Institute of Hydrology (BfG) and the Federal Waterways and Shipping Administration in Germany (WSA) in terms of existing monitoring stations or ships. All sampling campaigns were conducted between September 2020 and May 2022 at different discharge levels (high, moderate, and low discharge events).

Fig. 1 (A) Map of Europe with the location of the research area (red box). (B) Overview map of the Rhine basin with the different sampling locations in this study: the red circles indicate the locations of the fixed monitoring stations for the one-year data collection (Wesel/Rhine, Koblenz/Rhine, Koblenz/Moselle, and Weil/Rhine); the yellow triangles depict the locations of the sampling campaigns conducted from a ship or a pier (Koblenz (KAA), Niederlahnstein (LST), Brohl (BRL) and Emmerich (EMM)). The background map is a 90 m resolution Digital Elevation Model from the Shuttle Radar Topography Mission (SRTM GL3)³⁵ with the calculated stream network using MATLAB TopoToolbox by Schwanghart and Scherler (2014).³⁶

2.2 Sampling devices and materials used

To conduct river water sampling, in situ samplers were selected: 2 l Van-Dorn bottles (Aquatic BioTechnology, Spain), consisting of open tubes with openings and lids on both sides, were positioned along the flow direction of the river to minimize turbulence. Both lids close automatically by releasing a drop weight that hits a trigger. Moreover, a latex ring was ordered from the provider that is suitable for element analysis. Before their first utilization, they were leached in 2% HNO₃ (v/v) for 48 h, significantly improving the background concentrations of Cu, Zn and Sn.

For storage and transportation of the samples, 1 l PE sampling bottles (Kautex, Germany) were used. They were rinsed three times with 2% HNO₃ (v/v) and three times with ultrapure water (≤0.055 µS cm⁻¹ (corresponding to ≥ 18.2 MΩ cm)) before sampling. They were dried under clean room conditions (laminar flow box, SPETEC GmbH, Germany). PP vials (50 ml, VWR catalyst Laboratory Services, USA) were leached with 2% HNO₃ (v/v). After a minimum leaching time of 24 h, they were rinsed three times with ultrapure water and finally dried under clean room conditions.

For the filtration of the samples, 50 ml syringes (HENKE-JECT®, HENKE SASS WOLF GmbH, Germany) and syringe filters (0.45 µm, cellulose acetate, Minisart® NML, Sartorius, Germany) were used. All filters were pre-conditioned with the respective sample before filtration into 50 ml PP vials. The purification of all acids used for acidification of the filtered samples and for the digestion of the unfiltered samples – HNO₃ (65% w/w, EMSURE®, Merck GmbH, Germany) and HCl (37% w/w, EMSURE®, Merck GmbH, Germany) – was conducted via subboiling distillation (DST-1000, Savillex, USA).

2.3 Sampling campaigns

Four different sampling campaigns were conducted to cover a wide range of possible variabilities in river water monitoring with regard to multi-element analysis: (1) one-year temporal variability, (2) short-term temporal variability, (3) spatial variability and (4) small-scale spatial variability.

2.3.1 One-year temporal variability. All samples as well as discharge and turbidity data for the one-year temporal variability study originated from the environmental radioactivity monitoring network of the BfG. Automatic water samplers collect daily integrated river water samples over 24 h, creating 2.5 l samples which are stored at < 8 °C. The samples are collected and transported each month from the respective sampling location to the environmental radioactivity laboratory of the BfG. For the one-year data analysis, aliquots of the daily integrated samples from days 1 to 7 of every month from May 2021 to April 2022 were gathered. Not all daily integrated samples per month were assessed to reduce the total number of samples for this study. The following monitoring stations were considered for multi-element analysis (Fig. 1; n = 84 per sample type and station): Koblenz–Rhine (river-km: 590.3), Koblenz–Moselle (river-km: 2.0), Weil–Rhine (river-km: 170.3), and Wesel–Rhine (river-km: 814.0).

2.3.2 Short-term temporal variability. A high frequency sampling campaign was conducted at the stations Koblenz–Rhine (KAA), Brohl (BRL) and Emmerich (EMM) to assess possible short-term temporal variabilities (Fig. 1 and Table 1). Samples were taken at a fixed depth (25 cm from the surface) from a pier with two self-constructed 2 m sampling devices which contained a Van-Dorn bottle and a drop weight (Fig. S1). Samples were obtained at 1-minute intervals, to address the short-term temporal variability. All sampling campaigns were conducted between September 2020 and January 2022 (Table 1). Two high discharge campaigns were conducted at the station KAA, since the first sampling campaign (HW1: 5022 m³ s⁻¹) could not be performed from the pier due to inaccessibility. It had to be conducted from a pedestrian way with a railing for safety reasons. Therefore, a second high discharge campaign (HW2: 3528 m³ s⁻¹) was conducted from the pier at a lower discharge compared to HW1.

Table 1 Information on all sampling campaigns conducted for the short-term temporal, spatial and small-scale spatial variability in Koblenz (KAA), Niederlahnstein (LST), Brohl (BRL) and Emmerich (EMM) at the river Rhine. It has to be noted that the sampling station KAA was only used for the short-term temporal variability campaigns, while LST which is close to KAA was used for both spatial variability campaigns due to the ship availability at that location. In some cases and for different reasons (e.g., COVID pandemic, working safety restrictions, and ships not available), sampling was not feasible and the respective campaigns are marked with a “—”

		LST/KAA	BRL	EMM
	River-km (Rhine)	586.0/590.3	620.3	852.6
Short-term temporal variability		Low water level (LW)
	Sampling date	23.09.2020	23.09.2020	24.09.2020
	Discharge in m^3 s^−1	793	1130	1160
	No. of samples	11	10	11
		Moderate water level (MW)
	Sampling date	21.01.2021	21.01.2021	22.01.2021
	Discharge in m^3 s^−1	1417	1870	1590
	No. of samples	10	11	11
		High water level (HW)
	Sampling date	HW1: 02.02.2021	12.02.2021	10.01.2022
		HW2: 11.02.2021
	Discharge in m^3 s^−1	HW1: 5022	4560	4610
		HW2: 3528
	No. of samples	11/10	10	10
Spatial variability		Low water level (LW)
	Sampling date	17.11.2021	25.11.2021	05.10.2021
	Discharge in m^3 s^−1	811	850	1170
	No. of samples	56	52	66
		Moderate water level (MW)
	Sampling date	15.04.2021	22.04.2021	10.08.2021
	Discharge in m^3 s^−1	1255	1250	3160
	No. of samples	67	91	101
		High water level (HW)
	Sampling date	10.02.2022	—	10.01.2022
	Discharge in m^3 s^−1	2318	—	4610
	No. of samples	56	—	13
Small-scale spatial variability		Low water level (LW)
	Sampling date	30.11.2021	—	06.10.2021
	Discharge in m^3 s^−1	781	—	1170
	No. of samples	68	—	74
		Moderate water level (MW)
	Sampling date	—	29.11.2021	—
	Discharge in m^3 s^−1	—	1320	—
	No. of samples	—	76	—

---

2.3.3 Spatial variability. To assess spatial variabilities, cross-profile samplings were conducted with a self-constructed measuring frame that can hold five vertically fixed 2 l Van-Dorn bottles at 50 cm distances, thus covering a sampling depth of 2 m (Fig. S2). Sampling with the measuring frame was conducted from sampling ships provided by the Federal Waterways and Shipping Administration in Germany (WSA). These ships can target and remain fixed at precise sampling locations by anchoring themselves to the riverbed with hydraulic stamps and they have large cranes that allow sampling with heavy equipment. If these sampling ships were not available, e.g. at the sampling station in Emmerich, the ship was kept in place via GPS coordinates. Two 100 kg weights were additionally fixed to the bottom of the measuring frame to reduce the drifting of the frame within the river water current. The measuring frame was also fixed with a rope from the ship to avoid drifting or rotating around its own axis. To activate the automatic closing mechanism of the Van-Dorn bottles, a drop weight frame was constructed around the frame. The frame pressed on the closing mechanisms of the Van-Dorn bottles while released. Prior to sampling, the drop weight frame was lifted and fixed with a pin on a rope. Once the measuring frame reached the desired sampling location, the closing mechanism was triggered by pulling the rope and the drop weight frame caused the Van-Dorn bottles to close. The water depth during sampling was assessed through two different depth sensors, the YSI EXO 3 (Xylem Inc., USA) and the Nortek Vector (Nortek AS, Norway), which were attached to the measuring frame. Depending on the water depth, the measuring frame had to be deployed one to three times per cross-profile location to cover the whole river depth.

The spatial variability sampling campaigns were performed between April 2021 and February 2022 at the stations KAA, BRL and EMM (Table 1 and Fig. 1). Due to the limited availability of ships, especially during high discharge, not all discharge levels were covered in the project time. Two sampling campaigns at high discharge had to be aborted due to severe weather conditions and the associated safety and health risks. Moreover, due to the high discharge level during the sampling campaign at station LST (LST_HW), collection of a complete depth profile was not possible since the vertical measuring frame (with the additional weights) was not capable of withstanding the flow velocities at depth and the potential loss of the installation was foreseeable.

2.3.4 Small-scale spatial variability. In addition to the first described measuring frame, which will be referred to as the vertical measuring frame in the following, a second measuring frame with 16 2 l Van-Dorn bottles was constructed. The bottles were fixed in a 4 × 4 grid covering a 1.6 × 1.6 m window (Fig. S2). The 4 × 4 measuring frame was used to assess the small-scale spatial variability in the river cross-profile. The small-scale sampling was conducted at a defined depth from one riverbank to the other, but it did not cover the complete depth-profile since for this sampling campaign type, only the small-scale spatial variability within one grid is relevant.

2.4 Sample handling and storage

After sampling for the short-term temporal (Section 2.3.2), the spatial (Section 2.3.3), and small-scale spatial (Section 2.3.4) variability campaigns, an aliquot of 50 ml river water was filtered immediately using 50 ml syringes and syringe filters, which had been pre-conditioned with the respective sample. They were further acidified with purified HNO₃ to 2% HNO₃ (v/v). All filtered and acidified samples were stored frozen at −20 °C.

The unfiltered samples were stored at +4 °C in the dark in 1 l PE bottles after sampling. Samples from EMM were stored frozen (after filtration and acidification) or cooled (whole water samples) after transport to the BfG laboratory in Koblenz. Due to the large volume of 1 l whole water samples and capacity limitations, freezer storage was not feasible.

The one-year temporal variability samples were treated differently because they were derived from the regular German environmental radioactivity monitoring system of the BfG.³⁷ Within 1 to 7 days of arrival of the daily composite samples at the BfG laboratory, aliquots of 50 ml unfiltered water samples using 50 ml syringes and filtered aliquots of 50 ml (<0.45 µm) were stored frozen at −20 °C in 50 ml PP vials. Due to the unavoidable time delay in transport from the regular nationwide German monitoring network and as published previously,¹⁷ a potential influence of storage artifacts cannot be completely excluded for these samples.³⁸ However, this would most likely cause a systematic bias which does not affect the comparison of temporal variabilities within the one-year temporal dataset.

Before further analysis or sample preparation, all samples were acclimatized to room temperature.

2.5 Analytical protocols

All analyses were performed using an Agilent 8900 inductively coupled plasma triple quadrupole mass spectrometer ICP-QQQ-MS (Agilent Technologies, Japan) connected to an ESI autosampler (ESI SC4 DX, ESI Elemental Service & Instruments GmbH, Germany). Before each analysis, the instrument was tuned using 10 µg l⁻¹ of Li, Co, Y, Ce and Tl by checking the oxide formation rate (¹⁴⁰Ce¹⁶O⁺/¹⁴⁰Ce⁺ <2%), the doubly charged ion ratio (¹⁴⁰Ce⁺⁺/¹⁴⁰Ce⁺ <3%), the sensitivity, and the relative standard deviation (<3%) for all analysis modes/gases used (He, H₂, and O₂). Multi-element analyses were conducted according to Belkouteb et al.¹⁶ for filtered river water samples and Belkouteb et al.¹⁷ for unfiltered (whole water) river water samples. Both methods include 67 elements: Li, Be, B, (C), Na, Mg, Al, Si, P, S, (Cl), K, Ca, Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, Zn, Ga, Ge, As, Se, Br, Rb, Sr, Y, [Zr], Nb, Mo, Ru, Ag, Cd, In, Sn, Sb, Te, Cs, Ba, La, Ce, Pr, Nd, Sm, Eu, Gd, Tb, Dy, Ho, Er, Tm, Yb, Lu, Hf, Ta, W, Ir, Pt, [Au], Hg, Tl, Pb, Bi, Th, and U. For the filtered samples, Zr and Au, marked in brackets, were excluded from the method due to their proven instability in a HNO₃ matrix during method development and validation.¹⁶ However, Zr and Au are sufficiently stabilized for the analysis of whole water samples, while C and Cl cannot be analyzed due to the loss of C as CO₂ in the digestion process as well as the use of HCl as a digestion reagent.¹⁷ The whole water samples were digested using 8 ml 65% HNO₃ (w/w, sub-boiled) and 4 ml 37% HCl (w/w, sub-boiled) in a closed-vessel microwave system (MARS 6, CEM, Germany) with a ramp time of 20 min to 175 °C and a hold time of 20 min. For quality control of the digestion process, the sediment reference materials SdAR-M2 (International Association of Geoanalysts, UK) and NWWQB-1 (Environment and Climate Change Canada) were suspended in ultrapure water to simulate whole water samples as described by Belkouteb et al.¹⁷ Three quality control replicates as well as three blanks (ultrapure water) were digested in each digestion run (n = 40 including up to 31 samples per run). After digestion, the samples were completely transferred to volumetric flasks (50 ml) using a diluted digestion matrix solution (concentrated acid mixture of 2HNO₃ : 1HCl diluted in a 1 : 50 ratio with ultrapure water). They were subsequently analyzed for the element concentration via ICP-QQQ-MS. Each element was quantified in both sample types using a 10-point calibration. Further quality control measures comprise the direct analysis of surface water reference materials SPS-SW-1 and -2 (Spectrapure Standards AS, Norway) as well as NW-TMDA-51.5 (Environment and Climate Change Canada) by ICP-QQQ-MS and the use of the internal standards (ISTD) ¹⁰³Rh and ¹⁸⁵Re. For the following elements no certified reference material was available for the filtered surface water matrix: C, Cl, Ge, Br, Zr, Nb, Ru, In, Te, I, Hf, Ta, W, Ir, Pt, Au, and Hg. Therefore, calibration check solutions were used additionally for quality control. Details of the applied methods can be found in the respective publications.^16,17 After analysis and correction with the ISTD, the respective data sheet per measurement with all concentration data was exported from Agilent MassHunter Software (MassHunter 4.5, Workingstation, Version C01.05).

2.6 Statistical analysis

For further data analysis, R software (version 4.2.1, 2022 The R Foundation for Statistical Computing) was used. The limit of quantification (LOQ) was determined by using a minimum of ten blank replicates using either a 2% HNO₃ (v/v) solution for the filtered samples or a diluted digestion matrix solution (1 : 50) for the whole water samples. The following calculation was used: LOQ = x + 10 × SD with SD being the standard deviation and x being the mean of all blank concentrations per element. The results below the LOQ or with a relative standard deviation (RSD) above 10% for three analytical replicates per measurement were not considered in the further data analysis. Recoveries in a range of 80 to 120% were chosen as acceptable per element. The selection criteria for analytes from the total of 67 elements for detailed data analysis were as follows: (i) the availability of at least 50% of all considered data points or (ii) a visible trend at a certain event or campaign. Concentration–discharge-relationships (c–Q relationships) for the one-year data were determined via loglinear regression (base 10). Outliers in a row per element were detected based on the R package “outliers”, applying Grubbs' test with a statistical significance of p < 0.05. Pearson's correlation coefficients were calculated to determine correlations with p < 0.05 between elements and discharge or SPM using the R package “corrr” and “corrplot”. Standardization was performed with the R package “dplyr” by transforming the concentrations of one group into standardized values (z-scores) with a mean of 0 and a standard deviation (SD) of 1. K-means clustering using the elbow method was performed with the R packages “factoextra” and “stats” based on standardized concentration data (z-scores). By applying the elbow method to 69 samples, k = 2 was determined as the best number of clusters out of the available data.

3 Results and discussion

3.1 Temporal variability campaigns

3.1.1 One-year data. Additional information is available in the supporting information (SI) and at https://zenodo.org/records/17080553. The one-year dataset of the station Koblenz–Rhine includes the measured dissolved fraction (<0.45 µm), the measured whole water (unfiltered sample) and the calculated particulate concentration (difference between the two measured concentrations). K-means cluster analysis (elbow method; k = 2, n = 69) was conducted as an example on all whole water samples from the station Koblenz–Rhine to define element groups based on their similar behavior in varying element concentrations for one year. A clear pattern between two distinct element groups was observed for the whole water fraction (Fig. 2A): cluster 1 (P, V, Mn, Fe, Co, Ni, Cu, Ga, As, Rb, Y, Cs, Ba, La, Ce, Pr, Nd, Sm, Gd, Pb, and Th) and cluster 2 (B, S, K, Ca, Br, Sr, Mo, and U). Furthermore, a strong positive correlation with a Pearson correlation factor above 0.70 (p < 0.05) was found between discharge (Q) and V (0.80), Fe (0.81), Co (0.77), Ni (0.78), Ga (0.81), Rb (0.74), Y (0.77), Ce (0.74), Pr (0.73), Nd (0.74), Sm (0.71), Gd (0.77), and Th (0.71) in the whole water samples (Fig. 2B and Table S1). Conversely, for elements S (−0.80) and Br (−0.78), a strong negative correlation below −0.70 (p < 0.05) was identified. Therefore, in comparison to the cluster and correlation diagram, increasing discharge might lead to higher concentrations of elements in cluster 1 while the opposite is expected for cluster 2 (Fig. 2). The exemplary c–Q-relationships (Fig. 3) reflect the above-mentioned division of two element groups for Koblenz–Rhine within the whole water fraction. The higher the discharge, the higher the concentration for the elements Fe, Cu, Ga, As and Gd in the whole water fraction which are also part of cluster 1 (Fig. 2 and 3). Thus, a co-rising effect is observable. The opposite trend is observed for elements B and S which are grouped in cluster 2; hence, a dilution effect is observable at higher discharge in the whole water samples. At higher discharge, higher turbidity was measured due to higher SPM concentrations (Fig. S3). Consequently, elements showing a co-rising effect tend to be transported by SPM while the dilution effect elements are found predominantly in the dissolved fraction.

Fig. 2 Results of the one-year temporal variability campaign at Koblenz–Rhine of the whole water samples (first seven days of each month between May 2021 and April 2022, n = 84). (A) K-means clustering (k = 2, elbow method) of standardized data (z-scores of all one-year data per element) and (B) Pearson correlation matrix of selected elements and the discharge (Q). Exact r- as well as the respective p-values of the Pearson correlation matrix are listed in the SI (Tables S1 and S2). The colour indicates the sign of the r-value (−1 to 1) and the size of the dots corresponds to the value of the correlation coefficient (−1 to 0 and 0 to 1).

Fig. 3 Concentration–discharge (c–Q) plots (logarithmic scale, base 10) as log-linear regressions with a confidence interval of 95% (in a lighter shade) for the elements B, S, K, Fe, Ga, As, Ba, Gd as examples from the one-year dataset for water samples from the Koblenz–Rhine station (first seven days of each month between May 2021 and April 2022, n = 84). The particulate fraction was calculated as difference based on the available data for the filtered and whole water (unfiltered) fraction.

Exceptions to the above-mentioned pattern are K and Ca: almost constant regression lines are observed for the whole water fraction (Fig. 3 and S4) indicating a chemostatic behavior, which means that a variation in element concentration is absent or only slightly present in relation to the discharge.^19,31 Potassium (K) shows an almost zero correlation with discharge (r = 0.005, p < 0.05) confirming a chemostatic behavior in the whole water samples. Calcium (Ca) in the whole water samples shows highly significant positive correlations (p < 0.05) not only with elements grouped as dilution effect elements, e.g. with B (0.66), Sr (0.82) and U (0.89), but also with Ba (0.60) and a weak, but negative correlation with discharge (−0.20). Interestingly, even though Ba is grouped as a co-rising effect element in the cluster diagram, the correlation between discharge and Ba for the whole water samples is significantly positive (p < 0.05) but very low with r = 0.32 (Fig. 2). Significant positive correlations of Ba (with r above 0.6) with other elements that show a co-rising behavior, such as Fe (0.66) and Co (0.70) and the REEs, La (0.67), Ce (0.68), Pr (0.66), Nd (0.71), Sm (0.64), and Gd (0.69), as well as with elements grouped as dilution effect elements, such as K (0.79) and U (0.68), indicate multiple impacts and hence, a variable behavior of Ba (Fig. 2 and Table S1).

When comparing the different water fractions, Gd, for instance, shows opposite trends between the whole water and filtered samples at higher discharge: its concentration decreases in the filtered samples while it increases in the whole water samples (Fig. 3). Consequently, whole water Gd shows a co-rising effect while dissolved Gd is diluted at higher discharge. The dominance of particle-bound Gd in the whole water fraction is clearly visible when looking at the calculated particulate load in Fig. 3. This leads to the assumption that under high discharge levels several diffuse sources for particle-bound Gd exist in the catchment. The same pattern is observed for the elements Co, Ni, As, Rb, and Cs (Fig. 3 and S5 SI). Elements such as Fe and Ga show an almost constant or rather slightly increasing or decreasing concentration in the dissolved fraction while their concentration increases rapidly with higher discharge in the whole water fraction and, therefore, in the particle-bound fraction (Fig. 3). This pattern is also observable for P, V, Mn, Y, La, Ce, Pr, Sm and Pb (SI Fig. S5). The elements S, Br, and Mo show a dilution behavior in all fractions with increasing discharge (Fig. 3 and S5 SI) and thus, they may be interesting as indicators for other substances with a clearly defined dilution. For the filtered samples of the one-year data, more results are available for the stations Koblenz–Moselle, Weil–Rhine and Wesel–Rhine in the SI (SI Section 3). No further data exist for the whole water samples at Koblenz–Moselle, Weil–Rhine and Wesel–Rhine.

In the following sections, we divide the elements mainly into the previously established two groups: the dilution effect and co-rising effect elements. However, it is important to keep in mind that elements can deviate from this grouping, which is caused by multiple influencing factors. One example is Ba, which is grouped as a co-rising element in the cluster analysis, but shows a similar behavior to dilution as well as co-rising elements in the correlation diagram for the whole water samples (Fig. 2). Knapp et al.¹⁹ described a time-dependent effect in the Erlenbach catchment in Switzerland between May and November in 2017 and 2018, respectively. They reported that some elements and solutes in the dissolved fraction such as chloride, potassium and nitrate show a long-term dilution behavior, while on the event-scale (e.g. precipitation) the c–Q behavior varies between a dilution and co-rising pattern due to different reaction processes including atmospheric deposition or biological processes.¹⁹ This must also be considered when looking at long-term data. Different events may show varying effects on elements, which provides an opportunity to detect anomalies, e.g. if the input source is known at a certain discharge event, through precipitation, snow melt or anthropogenic sources. These findings indicate that long-term data are needed to define a baseline for each element at each sampling location or monitoring station, encompassing several annual cycles and events. By making anomalies visible at specific environmental events, the source of the element concentrations in rivers can be better determined if different water fractions are monitored, i.e. either the dissolved fraction together with SPM or the dissolved fraction together with the whole water samples. Pokrovsky et al.²¹ analyzed the dissolved fraction (<0.45 µm) of the Taz River in Western Siberia between 2015 and 2020. The elements Be, V, Cr, Co, Ni, Cu, Cd, Pb and Nb showed a positive correlation with river discharge for the dissolved fraction, which is also apparent in our study for the whole water samples; however, for Co and Ni, a dilution behavior is observed for the dissolved fraction underlining the importance of analyzing different river water fractions (Fig. 2). Interestingly, a strong increase in Fe concentration also becomes visible in the filtered fraction at higher discharge levels, by plotting the concentration against time, which is not as apparent as in Fig. 3 (refer to Fig. S6). Pokrovsky et al.²¹ found that Fe shows a varying transport behavior that could be influenced by several factors such as its mobilization under reducing conditions in groundwater or by organic colloids in surface waters. Two further examples of European Arctic rivers for which long-term multi-element data for the filtered fraction (<0.45 µm) are available are the Severnaya Dvina³⁹ and Pechora River.⁴⁰ They underline the varying behavior of Fe in Arctic rivers. Both studies revealed that three elemental groups for the filtered fraction could be defined depending on discharge and seasonal variations. For the Pechora River, for instance, generally decreasing concentrations with higher discharge were found for the elements Li, B, F, Na, Mg, Cl, Ca, Rb, Sr, Mo, Ba, W and U while co-rising behavior was observed for Be, Al, Ti, V, Ni, Ga, Cs, Se, Nb, Y, Zr, Hf, and Th and the REEs.⁴⁰ Varying behavior was observed in the third group including P, N, K, Fe, Mn, Cu, Zn and Mo.⁴⁰ Please note, that the river studied here is significantly different from Arctic rivers. When comparing the different water fractions, we found again that Rb exhibits a dilution behavior in the dissolved fraction in accordance with the findings for the Pechora River; however, in the whole water fraction, a co-rising behavior is observed underlining again the importance of the analysis of different water fractions to determine the total riverine element behavior. Another example is Ba, which showed varying behavior in our dataset unlike in the Pechora River for which a dilution effect was observed.⁴⁰

In addition to our results for the river Rhine, van der Perk and Vilches³³ found that the legacy elements Cr, Cu, Zn, Cd, Hg and Pb showed a negative relation with discharge in SPM due to the possible introduction of uncontaminated SPM at high discharge levels leading to a dilution effect in the riverine SPM. They analyzed SPM monitoring samples at the German–Dutch border in Lobith, which were gathered bi-weekly between 2011 and 2016.³³ Moreover, Liu et al.²² found that the partitioning between the sediment and water phase averaged over long time scales between 2000 and 2020 in the Yangtze River shows opposite trends for Cr (decreasing water and increasing sediment concentration) and Pb (increasing water and decreasing sediment concentration). Thus, long time series are necessary to identify changes in element concentrations between sample fractions and to detect anomalies over time. Details about the potential multiple impact factors could be detected if a long-term multi-element fingerprint dataset were available at the respective sampling location. The application of multi-element analysis is therefore a valuable tool to generate a long-term overview and for the understanding of a wide range of elements. However, for long-term data it is also important to gather complementary parameters such as pH, salinity or the organic matter content in addition to discharge and SPM concentration so that different impact factors can be better determined and interlinked. The importance of different influencing factors is highly dependent on the characteristics of a water body (e.g. catchment area and properties, climatic conditions, etc.) and thus the sampling strategy and analyte set should be chosen accordingly. The one-year data highlighted the already established possibility of grouping elements according to their behavior with discharge into co-rising, dilution, chemostatic and varying behavior.^19,31 However, the cluster analysis for the unfiltered water samples revealed that clustering the elements without additional parameters identifies mostly two groups of elements, which mainly show either a dilution or a co-rising behavior when compared to discharge variations. A varying behavior was only observable for Ba or K by calculating Pearson correlation coefficients. Moreover, the analysis of both, filtered and unfiltered river water, and the log-linear regression showed that the c–Q-relationships can differ between the sample types, e.g. for Co, Ni, As, Rb, Cs, and Gd. Since most of the available long-term element data refer to the dissolved fraction, long-term monitoring efforts should include the analysis of whole water samples in the future to complete the picture of element behavior in the water cycle.

3.1.2 Short-term temporal variability. To examine potential short-term temporal variations in element concentrations, the relative standard deviations (RSDs) of all samples taken per sampling campaign were compared. The RSDs of all 1 min intervals of the Koblenz–Rhine campaign (KAA) reveal mostly a higher variability in the filtered fraction (Fig. 4A). The RSDs of the filtered low discharge samples (KAA_LW) were mostly below 20%, except for Si, Yb, and Pb. However, all further sampling campaigns showed RSDs higher than 20%. This indicates that the short-term temporal variability increases with higher discharge in the dissolved fraction. The REEs, Ce, Pr, Nd, Sm, Tb, Dy, Ho, Er, Yb, and Lu, display one of the highest RSDs in the KAA_MW campaign of the filtered fraction.

Fig. 4 Relative standard deviation (RSD) in % of measured element concentrations for the complete dataset of the short-term temporal variability campaign for the station KAA at different discharge levels (LW, MW, and HW) for (A) filtered samples and (B) unfiltered samples. Samples were taken at 1-minute intervals for a total of 10 min (n = 10 or n = 11, Table 1). The dashed lines mark RSDs of 10% and 20%, respectively. For the KAA station, two high discharge campaigns (HW1 and HW2) were conducted as explained in Section 2.3.2.

An opposite trend is shown in Fig. 4B (whole water samples) where only a few elements show RSDs higher than 20%, such as Zn in three out of four sampling campaigns which is not observable at the further sampling locations (EMM and BRL, Fig. S7 and S8) indicating a possible point-source of Zn close to the KAA station. In addition, the REEs, Ce, Pr, Nd, and Er, again show higher RSDs between 20 and 30% in the KAA_MW campaign. These differences between the filtered fraction and the whole water samples were also visible at the other sampling locations (EMM and BRL, Fig. S7 and S8). Elements with high RSDs in the filtered fraction, for instance almost all REEs, V, Cr, Fe, Zn, Ga, Y or Pb (Fig. 4), are grouped as co-rising elements in Section 3.1.1; therefore with higher discharge, higher concentrations are expected for the co-rising elements in the whole water fraction.

Based on these results, it can be stated that the short-term temporal variability and the representativeness of a single sample is not an issue for whole water samples because of the low RSDs. It might, however, be important for the dissolved fraction. One explanation for the comparably higher RSDs of the filtered fraction can be a methodological uncertainty induced by the filtration procedure. With higher discharge, the particulate load increases and therefore the filter retention may change due to the formation of filter cakes, even though the filter was pre-conditioned with the sample and replaced up to five times per sample during high particulate load filtration to minimize these effects and to avoid filter rupture due to the increasing pressure during the filtration. As reported in the literature, filtration with a lower pore size is more suitable to determine the so-called “truly dissolved fraction” (<1 kDa) and to differentiate this fraction from the colloidal fraction.⁴¹ With such a bias, a short-term temporal variability cannot be clearly determined in the filtered fraction. However, the filtration was performed in a careful manner as already explained. Another explanation could be that differences in the 0.45 µm fraction are more pronounced than in the whole water fraction for the co-rising elements due to colloids or adsorption/desorption processes in the filtered fraction causing higher RSDs. As stated by Gaillardet et al.,⁴² sampling and filtration of organic-rich river samples can also encounter a pass-through of the colloidal phase which has to be taken into account when analyzing filtered samples. Consequently, if this is the case, short-term temporal variability cannot be excluded. These assumptions warrant further examination by using different filter pore sizes for the same sampling strategy applied in this study for future studies.

3.2 Spatial variability campaigns

3.2.1 Cross-profile measurements.
3.2.1.1 Hydrological aspects. The results of the spatial variability campaigns are influenced by the discharge during the campaign, the sampling location along the river, the position of the sampling bottle in the cross-profile of the river and lastly, also by possible temporal variabilities because the duration of one sampling campaign was up to 9 hours. The overview in Fig. 5A reflects, in accordance with Section 3.1.1, the two element groups: dilution and co-rising effect elements. The standardization was conducted per sampling station (LST, BRL, and EMM) and per element as explained in Section 2.6 (z-scores with mean = 0 and SD = 1). In the whole water samples, S, Ca, and U as well as B, Mg and Mo (not shown) are found in lower concentrations at high discharge levels (Fig. 5A). Exceptions are K and Br, as they generally follow the pattern of a dilution behavior with lower concentrations in comparison with the campaigns at low and moderate discharge levels (LW and MW). However, at high discharge (LST_HW and EMM_HW), their concentration increases (Fig. 5A). Furthermore, a high variability was apparent for LST_HW. It has to be noted at this point, that no sampling campaign for BRL during a high discharge event was possible as no ship was available at that time. When summarizing all data from the spatial sampling campaigns, K, Br and U, which are grouped as dilution effect elements (Fig. 5), show neither a strong positive nor a negative correlation with SPM (K:0.33, Br:0.24, and U:-0.21; p < 0.05).

Fig. 5 (A) and (B) Standardized concentrations (z-scores with mean = 0 and SD = 1) per element and per station (LST and BRL; EMM) sorted according to Section 3.1.1. Boxplots show the median as the middle line and the interquartile range (50% of all data) within the boxplots. (C) Pearson correlation diagram for the whole water samples from all spatial variability campaigns. Exact r- as well as the respective p-values of the Pearson correlation matrix can be found in the SI (Tables S4 and S5). The color indicates the sign of the r-value (−1 to 1) and the size of the dots corresponds to the value of the correlation coefficient (−1 to 0 and 0 to 1). For EMM_HW only one depth-profile was sampled due to the weather and run-off conditions during sampling.

For Al, Si, P, Ti, V, Mn, Fe, Ni, Cu, Ga, As, Rb, Y, Cs, La, Ce, Gd, Pb and Th, higher concentrations were determined at high discharge events with factors between 4 and 80, in comparison to a moderate or low discharge level at the LST and EMM stations. All elements in Fig. 5B clearly show a co-rising effect at the LST and EMM stations, except for Cu, Ba and Gd in EMM. Cu concentrations are highly variable at EMM_MW with the highest mean concentration of 12.9 µg l⁻¹ in comparison to EMM_LW (2.79 µg l⁻¹) and EMM_HW (6.04 µg l⁻¹) while Ba exhibits a dilution behavior at EMM and a co-rising effect at LST in the whole water fraction (Fig. 5B). Cu shows a Pearson correlation coefficient of only 0.06 with SPM (p < 0.05). The only element that shows a strong positive correlation with Cu is Si (0.85, p < 0.05; Fig. 5C and Table S4). According to van der Perk and Vilches (2020),³³ Cu along with other legacy elements, such as Cr, Zn, Cd, Hg, and Pb, originates from anthropogenic point-sources and these elements are mainly adsorbed on suspended particulate matter. The opposite behavior observed in our dataset for Koblenz–Rhine indicates that the contamination is not as pronounced as in other areas along the Rhine while, for instance, at the EMM station possible anthropogenic point-sources might be present (Fig. 5B).

Ba has a positive correlation with SPM (0.61, p < 0.05) and with elements which were previously grouped as co-rising elements (Al, Si, P, Ti, V, Mn, Fe, Ni, Ga, As, Rb, Y, Cs, La, Ce, Gd, and Pb: r > 0.6, p < 0.05) as well as with the elements K and U (r > 0.6, p < 0.05; Fig. 5C and Table S4), which indicates a varying behavior of Ba in the whole water samples. Interestingly, Gd generally shows a higher concentration in EMM than at all other sampling locations, except for the campaign LST_HW with the highest overall concentration, suggesting a distinct point source upstream of EMM (Fig. 5B). Moreover, for the BRL and LST stations a decrease in element concentration in the whole water samples was observed between LW and MW (Fig. 5B), which could be due to the minimal differences in discharge of the two campaigns (refer to Table 1).

When looking at both the one-year temporal and the spatial data, it is evident that K and Ba exhibit highly variable concentrations in both datasets. While K shows no correlation with discharge (0.005, p < 0.05) or a minimal positive correlation with SPM (0.33, p < 0.05), it shows the highest correlation with Ba in both datasets (temporal: 0.79; spatial: 0.74; p < 0.05). For K, variable behavior was also observed in other surface waters, such as the Pechora River⁴⁰ or in the Erlenbach catchment.¹⁹ K is generally a mobile element. However, in the streams studied in the Erlenbach catchment, for example, it showed a long-term dilution behavior, but weathering processes potentially generate particulate-bound K which follows a co-rising behavior. Ba shows a significant positive correlation with SPM (0.61, p < 0.05) and a minimal positive correlation with discharge (0.32, p < 0.05). Coffey et al.⁴³ reported that the fractionation change from particulate Ba to the dissolved fraction is highly dependent on the SPM concentration, the salinity and the river discharge. At high salinity levels, and consequently, higher concentrations of Na or K, and low discharge, more dissolved Ba is available than at high discharge. In laboratory experiments, a release from the particulate to the dissolved fraction occurred at low salinity. This could potentially explain some aspects of the highly variable behavior of Ba in this study, although studies from freshwater environments and marine waters can only be compared to a limited extent. Moreover, Coffey et al.⁴³ also mentioned a high dependency of Ba on Mn and Fe due to possible co-precipitation with their (oxyhydr)oxides. Positive correlations were observed in this study for Ba with Mn and Fe above 0.6 (p < 0.05) in the whole water fraction (refer to Fig. 2 and 5), which emphasizes a possible interlinked behavior with the co-rising effect elements Fe and Mn.

The exceptional element behavior in this dataset of K, Br, and U as well as Cu, Gd and Ba underlines again that different influencing factors are relevant for each sampling location, time, and discharge.

3.2.1.2 Spatial aspects. The application of the vertical measuring frame gathering five samples at a time at a defined depth revealed lateral variabilities from one riverbank to the other within the river cross-section (Fig. 6 and 7). The elements Gd, K, P and Fe, shown as examples, cover different factors and element groups: Gd is an REE and exhibited a distinct behavior at EMM due to a possible point-source, as indicated in the previous section, K is an alkali metal and showed variable behavior in the previous datasets, P is a non-metal, and Fe is a commonly analyzed element in river monitoring. Further results are shown in SI Section 6.

Fig. 6 Concentrations in µg l⁻¹ of Gd and K in the spatial variability sampling campaigns EMM_MW (Emmerich–Rhine at a moderate discharge level (A–C and G–I)) and EMM_LW (Emmerich–Rhine at a low discharge level (D–F and J–L)) in different sampling fractions: filtered fraction (top row), whole water samples (middle row) and calculated particulate load (bottom row). White circles indicate measurements below the LOQ or with an analytical RSD > 10%.

Fig. 7 Concentrations in µg l⁻¹ of P and Fe in the spatial variability sampling campaigns EMM_MW (Emmerich–Rhine at a moderate discharge level (A–C and G–I)) and EMM_LW (Emmerich–Rhine at a low discharge level (D–F and J–L)) in different sampling fractions: filtered fraction (top row), whole water (middle row) and calculated particulate load (bottom row). White circles indicate measurements below the LOQ or with an analytical RSD > 10%.

Detailed investigations for Gd, as one example, show that two factors play a major role when higher concentrations are observed in a river cross-section from the left to the right river bank, e.g. in EMM_MW with a factor of around 2 (Fig. 6): (i) a possible anthropogenic point-source upstream and close to EMM, leading to higher concentrations in the dissolved fraction and (ii) the association of REEs with SPM in the whole water samples with a strong positive correlation of La, Ce and Gd with SPM for all spatial variability campaign data (0.93, 0.96, 0.83, p < 0.05; Fig. 5C and SI Table S4). The strong association of REEs with SPM was also observed for the Moselle by Le Meur et al.⁴⁴ Furthermore, Gd is an element for which an anomaly, in comparison to other REEs, has already been described not only in the Rhine but also worldwide in industrialized areas.^10,13 Kulaksız and Bau¹³ detected not only Gd anomalies in the Rhine, but also anomalies for La and Sm, caused by anthropogenic point sources. In the case of Gd, these anomalies are traced back to Gd compounds used in MRI (magnet resonance imaging) contrast agents.^13,14 It is known that Gd-based contrast agents are only transported in the dissolved fraction,^13,45 which could be the case for the sampling campaign in Emmerich at a moderate water level (EMM_MW, Fig. 6A). In comparison to the sampling campaign in Emmerich at a low water level (EMM_LW, Fig. 6D), Gd concentrations are around 300 times higher in the dissolved fraction for EMM_MW. As observed in Section 3.1.1, a higher Gd concentration in the dissolved fraction at a lower discharge would have been expected in the long-term data. However, there is no continuous monitoring of Gd at the sampling station in Emmerich. Based on the spatial variability campaign data, it has to be assumed that an introduction of dissolved Gd species into the river occurred at the sampling time of EMM_MW (moderate water level) but not of EMM_LW (low water level). Regarding the whole water samples, higher concentrations of SPM also lead to higher concentrations of particulate bound Gd (Fig. 5C). For EMM_LW, an increase in SPM concentrations from the left to the right river side up to a factor of 4 was observed, from 5.00–10.9 mg l⁻¹ at the left river side to 16.3–23.8 mg l⁻¹ at the right river side, which explains the higher Gd concentration in the whole water fraction compared to the dissolved fraction. For K, an increase in the particulate fraction is evident in both campaigns EMM_MW and EMM_LW from the left to the right river side, as shown exemplarily in Fig. 6I and L.¹⁹ The opposite pattern was observed in the filtered fraction (Fig. 6G and J) which underlines the importance of the analysis of several sample fractions to gain more detailed insights. P shows the same pattern as K for EMM_MW but not for EMM_LW (Fig. 7). At low water levels, P concentrations are evenly distributed across all fractions in the river cross-profile. Fe seems to show no spatial pattern for the whole water fraction having a concentration 100 times higher than in the filtered fraction for EMM_MW. However, in the filtered fraction a clear decrease in concentration from the left to the right riverbank of a factor 2 was observed for EMM_MW and EMM_LW (Fig. 7G and J).

Depth-profile variations in element concentrations were less apparent in our datasets than lateral differences within the river cross-profile. Spatial variabilities are known to be more pronounced with depth than across the river.^24,26 Lupker et al. (2011)²⁶ showed for the Ganga basin, that while Na and Si concentrations decreased in SPM, Al, Fe, Mg, and, K concentrations increased with depth because of different settling velocities. The examination of SPM monitoring data within our project (20+ years) revealed similar results for SPM concentrations, indicating that lateral variabilities are more pronounced than vertical variabilities.⁴⁶ According to Slabon et al.,⁴⁶ an increased lateral variability in the Rhine is more likely due to turbulent mixing processes and, to a lesser extent, affected by site specific properties (e.g. channel curvature and width or sedimentary characteristics). In contrast to our study, Bouchez et al.²⁴ found that an underestimation might occur when samples are only grabbed from the river surface by conducting cross-section measurements at the Amazon River, e.g. for Al and Cs in the particulate fraction by a factor 2. They attributed this to the hydrodynamic sorting of larger grain sizes towards the bottom. Even though grain size measurements were not included in this study, it is evident that the lateral variability is stronger than the depth-profile variability at the Rhine. This could be explained by the fact that the Rhine is rather shallow (e.g., the highest depth during the sampling campaigns was around 7 m) compared to the Amazon River (depths up to 60 m,).²³ In shallow rivers, different mechanisms such as river confluences might drive a stronger lateral variability than a depth-profile variability. Moreover, the water flow might be less turbulent resulting in a stronger lateral gradient.

The insights obtained from the spatial sampling campaigns not only highlight the variability between different discharge levels, but also the variable element distribution in river cross-profiles. This is an important aspect when it comes to possible uncertainties in monitoring data, especially since monitoring data refer mostly to a single measurement point, especially in regulatory approaches. Unlike further publications on spatiotemporal variabilities in the river cross-profile, uncertainties in river monitoring of the Rhine are more pronounced on the lateral scale than with depth. Differences in element concentrations with factors of up to 4 between the right and left riverbank indicate that single-point sampling cannot provide a comprehensive picture regarding element behavior in rivers.

3.2.1.3 Spatiotemporal aspects. Fig. 6 and 7 indicate again that short-term temporal variabilities could play a role. Overlapping sampling points are available because the vertical measuring frame only covers a depth of 2 m with 5 simultaneous samples, and hence, it had to be positioned several times in the depth-profile to cover the whole profile from the river surface to the riverbed. During this sampling procedure, some of the sampling bottles overlapped with samples taken previously. The overlapping points in Fig. 6 and 7, e.g. for K or Fe in EMM_LW, show concentration differences in the filtered and whole water fractions corresponding to about a 26% change for K in the filtered fraction or to about 37% for Fe in the whole water fraction, respectively. This is a good example of possible short-term temporal variations. The river continues flowing and a different time stamp is recorded for the second sampling point. In Fig. 4 and Section 3.1.2., it is shown that short-term temporal variations were rather low when analyzing several samples (n = 10) taken at 1-minute intervals from the same spot for the whole water samples, but relatively high in the filtered fraction suggesting short-term temporal variability in element concentration. Consequently, spatial variabilities could be also an expression of temporal variabilities.

3.2.2 Small-scale spatial variability. The small-scale spatial variability covers distances of 40 cm between each sampling point in a 4 × 4 grid with a total vertical and lateral length of 1.6 m × 1.6 m. This frame was positioned 5 times in the river from the left to the right riverbank (L01 to L05) at intermediate depth. The results from one sampling campaign in Emmerich during low water level are presented as one example in Fig. 8. Further results are depicted in SI Section 7.

Fig. 8 Concentrations of Gd and Ga in µg l⁻¹ in the small-scale spatial variability sampling campaigns in EMM in different sampling fractions: filtered fraction (A and D), whole water (B and E) and calculated particulate fraction (C and F). L01 to L05 indicate one 4 × 4 sampling grid (16 sampling bottles) at the different sampling locations in the river cross-profile with the following distances from the left riverbank: L01: 35.5 m, L02: 105.5 m, L03: 175.5 m, L04: 245.5 m, and L05: 315.5 m. Grey spots indicate measurements below the LOQ or with RSDs > 10%. White spots indicate sampling spots which could not be sampled.

Small-scale variabilities within one sampling grid were not as pronounced as the cross-section differences for the station EMM. They were observed, for example, for L05 in the calculated particulate fraction for Gd (Fig. 8C) with a factor of around 1.7 between the maximum concentration (149 ng l⁻¹) to the mean concentration (85.8 ng l⁻¹). For all further element concentrations in the whole water and filtered fractions, variations within a sampling grid were below a factor of 2 for all elements (refer to SI Section 7).

Both low discharge sampling campaigns of the vertical and small-scale measuring frame in Emmerich were conducted at subsequent days and similar discharge (refer to Table 1). It is again emphasized that the Gd concentration is high in the dissolved fraction at the sampling station EMM with concentrations between 180 and 323 ng l⁻¹ at the right river side (L04 and L05 for Gd in Fig. 8) which is the dominant fraction in the whole water samples (265 to 389 ng l⁻¹) while the calculated particulate concentration is lower with concentrations between 16.7 and 174 ng l⁻¹. Several entry paths of anthropogenic Gd could be available at this sampling point, which is exclusively present in the dissolved fraction and most prominent in effluents of wastewater treatment plants (WWTPs) containing Gd-based MRI contrast agents.¹³ The sampling location at the right side of the Rhine in Emmerich is influenced by several metropolitan areas along the Rhine (Cologne, Düsseldorf and the major part of the Ruhr basin with about 5 million inhabitants) and the respective right tributaries Ruhr, Emscher and Lippe with distances of around 40 km, 50 km and 80 km, respectively, from the sampling location in Emmerich. The high concentration of dissolved Gd near the right riverbank indicates that a nearby point source could be present, or that the mixing of river confluences was not yet completed at the position of the sampling location.^46,47 Therefore, a tributary or a WWTP could only be detected at one river side due to insufficient mixing of the two different water masses.

Additionally, the rarely analyzed elements Ga (Fig. 8) and Y (refer to SI Section 7), both grouped as co-rising elements, show higher concentrations on the left river side (L01, L02 and L03) for the whole water samples (Fig. 8 and SI Section 7). Furthermore, the particulate fraction dominates the total Ga load with factors of 10 to 20 between the particulate and dissolved fractions. Overall, the small-scale spatial variability campaign underlines that lateral variabilities in element concentration within the river cross-section are more pronounced than small-scale spatial variabilities.

3.3 Overall challenges and lessons learned from the study

The complex and cross-disciplinary nature of this study was associated with several logistical and technical challenges referring to e.g., sampling, storage, and analysis. These were already mentioned in the respective section but are presented in a condensed manner here once again. (I) Time pressure arising from heavy weather events: the reaction time for high discharge events was very short and all safety requirements had to be considered as high discharge bears several occupational safety and health risks. For future studies, we advise that it is very important to maintain a high level of communication between scientists and ship crews, with clear decision structures.

(II) Many different samples were collected and analyzed in this study. Due to various logistical and organizational constraints, not all samples could be processed in exactly the same manner. However, these differences are clearly described and resulting limitations are openly discussed in the manuscript. For example, since the one-year temporal variability data were derived from the regular monitoring network of the BfG, time delay in processing the samples was unavoidable due to the long transport distances to the main laboratory.

(III) A general lack of certified multi-element surface water reference materials for several elements remains a challenge in the application of multi-element analyses, which was discussed already by Belkouteb et al.,^16,17 providing a suitable workaround by using in-house quality control checks.

(IV) When analyzing 67 elements, it is evident that not all elements are present above the respective LOQs in the respective river. Therefore, for further data analysis we defined that at least 50% of the data must be above the LOQ for the next data processing steps. In addition, as an additional quality measure, element concentrations were not considered if recoveries per element in the quality checks were outside of the acceptable range of 80 to 120%.

4 Conclusions

Multi-element analysis in a single analytical run enables an almost complete and fast overview of the element content in river water samples. Two methods yielding concentrations of 67 elements for filtered (<0.45 µm) and unfiltered samples (whole water samples) were used to gain a detailed picture of the element distribution on spatial and temporal scales at three locations in the German part of the Rhine (Koblenz, Brohl, and Emmerich) and different discharge levels (low, moderate and high). The temporal and spatial variability campaign data showed in relation to discharge that the elements can be mainly grouped into two main groups within the whole water fraction: elements showing a dilution effect (e.g.: B, Mg, S, K, Ca, Br, Sr, Mo, and U) or a co-rising effect (e.g.: Al, Si, P, Ti, V, Mn, Fe, Co, Ni, Cu, Ga, As, Rb, Y, Cs, Ba, La, Ce, Pr, Nd, Sm, Gd, Pb, and Th) with higher discharge. This grouping is generally used as a starting point to detect anomalies through variable element concentrations caused by specific anthropogenic and environmental effects at a certain spot, which can be beneficial e.g. for AI training datasets in the future. The investigation of element concentrations in both the unfiltered and the filtered samples over the one-year temporal (May 2021 to April 2022) or the spatial scale (e.g. Emmerich with around 350 m river width and 7 m river depth) in relation to discharge or SPM content, allowed the identification and description of further differences between elements and their partitioning, such as for K, Cu, Ba or Gd. The analysis of filtered as well as unfiltered river water samples allowed us to observe opposite long-term element behavior between the fraction, e.g. dissolved Gd showed a dilution effect while whole water Gd followed a co-rising effect which was also visible for Co, Ni, As, Rb, and Cs. This observation emphasizes the need for the analysis of different fractions, e.g. to analyze at least either the dissolved and whole water fraction or the dissolved and particulate fraction to provide a better understanding of the processes at monitoring sites.

The simultaneous sampling of five samples in a vertical row in 50 cm distances for assessing spatial variabilities while minimizing the effect of temporal variabilities revealed that depth gradients are not as relevant as lateral variabilities for the Rhine which is a rather shallow river compared to the large rivers for which depth gradients were previously described in cross-profile measurements.²³ In our study, lateral differences in element concentrations of up to a factor of 4 were observed which is likely due to delayed mixing of water masses after the inflow of tributaries or due to possible point-sources located close to the sampling points. The novel sampling strategy for small-scale variabilities developed for this study provided insight into possible short-term changes; however, these variations can be considered to be only of minor importance for the whole water fractions, while it might be interesting to follow up on more detailed investigations for the filtered fraction.

In summary, the following main conclusions can be drawn from our different types of sampling campaigns:

• Long-term multi-element data with different sample fractions (e.g., at least filtered and unfiltered water) are necessary to comprehensively detect anomalies over time.

• Short-term temporal variabilities were not observed in the whole water fraction while they can be present to some extent in the filtered fraction (<0.45 µm).

• Spatial variabilities in river cross profiles are dominated by lateral gradients (at least for shallow rivers such as the Rhine) while depth gradients were less pronounced.

• Cross-profile differences are more pronounced than small-scale spatial variabilities.

Consequently, sampling at only one river side may induce uncertainties in measured concentrations and single-point monitoring activities should include an initial and ongoing site-specific characterization to better describe potential blind spots. As an example, regular sampling of both river sides as a calibration tool could be an option to reduce uncertainties in monitoring, especially if catchment wide statements such as the establishment of element budgets are required. The presence of possible point-sources such as Gd or Cu indicates the necessity for frequent monitoring with sufficient spatial resolution to detect anthropogenic impacts, especially in densely populated areas.

Conflicts of interest

There are no conflicts to declare.

Data availability

All data are found in tables or figures within the manuscript, the supplementary information (SI) and at https://zenodo.org/records/17080553. If further information is required, it can be requested from the corresponding author. Supplementary information: additional information on methods and materials, graphs and tables on the one-year temporal data, the short -term temporal variability data, the spatial variability campaigns and the heat maps of the campaigns; the repository contains the data in tables. See DOI: https://doi.org/10.1039/d5em00487j.

Acknowledgements

This study was funded by the German Federal Ministry of Transport. The study was part of the research project “URSACHEN” (“Investigating the uncertainties of spatiotemporal variable substance load estimations in rivers”). We would like to thank all colleagues of the Federal Institute of Hydrology (BfG) who helped during the sampling campaigns and in laboratory work: Frank Ruenz, Rolf Hartmann, Ulf Koeppe, Tim Fahrenkrog, Pia Parrhysius, Ingo Habermann, Sascha Koch, Stefan Kiefer. Furthermore, we are very thankful for all the help provided by the different ship crews of the “Walter Türk” and “VSS Emmerich” of the WSA. We thank the anonymous reviewers for their helpful comments.

References

1. M. T. H. van Vliet, J. Thorslund, M. Strokal, N. Hofstra, M. Flörke, H. Ehalt Macedo, A. Nkwasa, T. Tang, S. S. Kaushal, R. Kumar, A. van Griensven, L. Bouwman and L. M. Mosley, Global river water quality under climate change and hydroclimatic extremes, Nat. Rev. Earth Environ., 2023, 4, 687–702.
2. E. R. Jones, D. J. Graham, A. van Griensven, M. Flörke and M. T. H. van Vliet, Blind spots in global water quality monitoring, Environ. Res. Lett., 2024, 19, 091001.
3. A. Najah Ahmed, F. Binti Othman, H. Abdulmohsin Afan, R. Khaleel Ibrahim, C. Ming Fai, M. Shabbir Hossain, M. Ehteram and A. Elshafie, Machine learning methods for better water quality prediction, J. Hydrol., 2019, 578, 124084.
4. W. Zhi, A. P. Appling, H. E. Golden, J. Podgorski and L. Li, Deep learning for water quality, Nat. Water, 2024, 2, 228–241.
5. Q. Duan, Q. Zhang, X. Quan, H. Zhang and L. Huang, Innovations of water pollution traceability technology with artificial intelligence, Earth Crit. Zone, 2024, 1, 100009.
6. S. Huang, Y. Wang and J. Xia, Which riverine water quality parameters can be predicted by meteorologically-driven deep learning?, Sci. Total Environ., 2024, 946, 174357.
7. G. E. Batley and P. G. C. Campbell, Metal contaminants of emerging concern in aquatic systems, Environ. Chem., 2022, 19, 23–40.
8. A.-M. Dendievel, C. Grosbois, S. Ayrault, O. Evrard, A. Coynel, M. Debret, T. Gardes, C. Euzen, L. Schmitt, F. Chabaux, T. Winiarski, M. Van Der Perk and B. Mourier, Key factors influencing metal concentrations in sediments along Western European Rivers: A long-term monitoring study (1945–2020), Sci. Total Environ., 2022, 805, 149778.
9. L. Duester, N. Livrozet, S. Poturalski, T. Stoetter, A.-L. Gerloff and M. D. Heintz, The ICPR measuring programme chemistry and its monitoring approach – a look back and a glimpse of the future, Water Int., 2024, 49, 446–454.
10. M. Bau and P. Dulski, Anthropogenic origin of positive gadolinium anomalies in river waters, Earth Planet. Sci. Lett., 1996, 143, 245–255.
11. S. Kulaksız and M. Bau, Contrasting behaviour of anthropogenic gadolinium and natural rare earth elements in estuaries and the gadolinium input into the North Sea, Earth Planet. Sci. Lett., 2007, 260, 361–371.
12. S. Kulaksız and M. Bau, Rare earth elements in the Rhine River, Germany: First case of anthropogenic lanthanum as a dissolved microcontaminant in the hydrosphere, Environ. Int., 2011, 37, 973–979.
13. S. Kulaksız and M. Bau, Anthropogenic dissolved and colloid/nanoparticle-bound samarium, lanthanum and gadolinium in the Rhine River and the impending destruction of the natural rare earth element distribution in rivers, Earth Planet. Sci. Lett., 2013, 362, 43–50.
14. P. Louis, D. A. L. Vignati, S. Pontvianne and M.-N. Pons, Spatial distribution of Rare Earth Elements in a transnational watershed: The case of the Danube River, Sci. Total Environ., 2023, 892, 164368.
15. M. Filella, I. Pomian-Srzednicki and P. M. Nirel, Development of a powerful approach for classification of surface waters by geochemical signature, Water Res., 2014, 50, 221–228.
16. N. Belkouteb, H. Schroeder, J. Arndt, J. G. Wiederhold, T. A. Ternes and L. Duester, Quantification of 68 elements in river water monitoring samples in single-run measurements, Chemosphere, 2023, 320, 138053.
17. N. Belkouteb, H. Schroeder, J. G. Wiederhold, T. A. Ternes and L. Duester, Multi-element analysis of unfiltered samples in river water monitoring—digestion and single-run analyses of 67 elements, Anal. Bioanal. Chem., 2024, 416, 3205–3222.
18. M. A. Gunn, Z. S. Moran, B. M. Pracheil and P. J. Allen, Spatial changes in trace element and strontium isotope water chemistry in a temperate river system with application to sturgeon movement, J. Freshw. Ecol., 2019, 34(1), 739–755.
19. J. L. A. Knapp, J. von Freyberg, B. Studer, L. Kiewiet and J. W. Kirchner, Concentration–discharge relationships vary among hydrological events, reflecting differences in event characteristics, Hydrol. Earth Syst. Sci., 2020, 24, 2561–2576.
20. R. R. Pant, F. Zhang, F. U. Rehman, M. Koirala, K. Rijal and R. Maskey, Spatiotemporal characterization of dissolved trace elements in the Gandaki River, Central Himalaya Nepal, J. Hazard. Mater., 2020, 389, 121913.
21. O. S. Pokrovsky, R. M. Manasypov, A. V. Chupakov and S. Kopysov, Element transport in the Taz River, western Siberia, Chem. Geol., 2022, 614, 121180.
22. Q. Liu, X. Xu, L. Lin, L. Bai, M. Yang, W. Wang, X. Wu and D. Wang, A retrospective analysis of heavy metals and multi elements in the Yangtze River Basin: Distribution characteristics, migration tendencies and ecological risk assessment, Water Res., 2024, 254, 121385.
23. J. Bouchez, J. Gaillardet, C. France-Lanord, L. Maurice and P. Dutra-Maia, Grain size control of river suspended sediment geochemistry: Clues from Amazon River depth profiles, Geochem., Geophys., Geosyst., 2011, 12, Q03008.
24. J. Bouchez, M. Lupker, J. Gaillardet, C. France-Lanord and L. Maurice, How important is it to integrate riverine suspended sediment chemical composition with depth? Clues from Amazon River depth-profiles, Geochim. Cosmochim. Acta, 2011, 75, 6955–6970.
25. J. Bouchez, V. Galy, R. G. Hilton, J. Gaillardet, P. Moreira-Turcq, M. A. Pérez, C. France-Lanord and L. Maurice, Source, transport and fluxes of Amazon River particulate organic carbon: Insights from river sediment depth-profiles, Geochim. Cosmochim. Acta, 2014, 133, 280–298.
26. M. Lupker, C. France-Lanord, J. Lavé, J. Bouchez, V. Galy, F. Métivier, J. Gaillardet, B. Lartiges and J.-L. Mugnier, A Rouse-based method to integrate the chemical composition of river sediments: Application to the Ganga basin, J. Geophys. Res.: Earth Surf., 2011, 116, F04012.
27. M. Ghotbizadeh, C. W. Cuss, I. Grant-Weaver, A. Markov, T. Noernberg, A. Ulrich and W. Shotyk, Spatiotemporal variations of total and dissolved trace elements and their distributions amongst major colloidal forms along and across the lower Athabasca River, J. Hydrol. Reg. Stud., 2022, 40, 101029.
28. V. V. Gordeev, O. S. Pokrovsky, A. V. Zhulidov, A. S. Filippov, T. Y. Gurtovaya, R. M. Holmes, L. S. Kosmenko, J. W. McClelland, B. J. Peterson and S. E. Tank, Dissolved Major and Trace Elements in the Largest Eurasian Arctic Rivers: Ob, Yenisey, Lena, and Kolyma, Water, 2024, 16, 316.
29. J. W. McClelland, R. M. Holmes, B. J. Peterson, R. Amon, T. Brabets, L. Cooper, J. Gibson, V. V. Gordeev, C. Guay, D. Milburn, R. Staples, P. A. Raymond, I. Shiklomanov, R. Striegl, A. Zhulidov, T. Gurtovaya and S. Zimov, Development of a Pan-Arctic Database for River Chemistry, Eos Trans. Am. Geophys. Union, 2008, 89, 217–218.
30. P. Floury, J. Gaillardet, E. Gayer, J. Bouchez, G. Tallec, P. Ansart, F. Koch, C. Gorge, A. Blanchouin and J.-L. Roubaty, The potamochemical symphony: new progress in the high-frequency acquisition of stream chemical data, Hydrol. Earth Syst. Sci., 2017, 21, 6153–6165.
31. M. Bieroza, S. Acharya, J. Benisch, R. N. ter Borg, L. Hallberg, C. Negri, A. Pruitt, M. Pucher, F. Saavedra, K. Staniszewska, S. G. M. van’t Veen, A. Vincent, C. Winter, N. B. Basu, H. P. Jarvie and J. W. Kirchner, Advances in Catchment Science, Hydrochemistry, and Aquatic Ecology Enabled by High-Frequency Water Quality Measurements, Environ. Sci. Technol., 2023, 57, 4701–4719.
32. M. Harhash, J. G. Wiederhold, I. Nett, J. Arndt, A. Zavarsky, T. A. Ternes and L. Duester, Automated high-temporal-resolution multi-element river water monitoring with ICP-MS, Environ. Sci. Eur., 2025, 37, 89.
33. M. van der Perk and A. E. Vilches, Compositional dynamics of suspended sediment in the Rhine River: sources and controls, J. Soils Sediments, 2020, 20, 1754–1770.
34. R. S. E. W. Leuven, G. van der Velde, I. Baijens, J. Snijders, C. van der Zwart, H. J. R. Lenders and A. bij de Vaate, The river Rhine: a global highway for dispersal of aquatic invasive species, Biol. Invasions, 2009, 11, 1989–2008.
35. NASA Shuttle Radar Topography Mission (SRTM), 2013.
36. W. Schwanghart and D. Scherler, Short Communication: TopoToolbox 2 – MATLAB-based software for topographic analysis and modeling in Earth surface sciences, Earth Surf. Dyn., 2014, 2, 1–7.
37. G. Dersch and W. J. Krause, Monitoring and prediction of the dispersion of radioactive substances in German Federal waterways – concepts and methods, Kerntechnik, 2004, 69, 214–222.
38. A. Przibilla, S. Iwainski, T. Zimmermann and D. Pröfrock, Impact of storage temperature and filtration method on dissolved trace metal concentrations in coastal water samples, Water Environ. Res., 2023, 95, e10922.
39. A. V. Chupakov, O. S. Pokrovsky, O. Y. Moreva, L. S. Shirokova, N. V. Neverova, A. A. Chupakova, E. I. Kotova and T. Y. Vorobyeva, High resolution multi-annual riverine fluxes of organic carbon, nutrient and trace element from the largest European Arctic river, Severnaya Dvina, Chem. Geol., 2020, 538, 119491.
40. A. V. Chupakov, O. S. Pokrovsky, O. Y. Moreva, E. I. Kotova, T. Y. Vorobyeva and L. S. Shirokova, Export of organic carbon, nutrients and metals by the mid-sized Pechora River to the Arctic Ocean, Chem. Geol., 2023, 632, 121524.
41. C. Guéguen and J. Dominik, Partitioning of trace metals between particulate, colloidal and truly dissolved fractions in a polluted river: the Upper Vistula River (Poland), Appl. Geochem., 2003, 18, 457–470.
42. J. Gaillardet, J. Viers and B. Dupré, in Treatise on Geochemistry, Elsevier, 2014, pp. 195–235.
43. M. Coffey, F. Dehairs, O. Collette, G. Luther, T. Church and T. Jickells, The Behaviour of Dissolved Barium in Estuaries, Estuarine, Coastal Shelf Sci., 1997, 45, 113–121.
44. M. Le Meur, E. Montargès-Pelletier, A. Bauer, R. Gley, S. Migot, O. Barres, C. Delus and F. Villiéras, Characterization of suspended particulate matter in the Moselle River (Lorraine, France): evolution along the course of the river and in different hydrologic regimes, J. Soils Sediments, 2016, 16, 1625–1642.
45. G. Klaver, M. Verheul, I. Bakker, E. Petelet-Giraud and P. Négrel, Anthropogenic Rare Earth Element in rivers: Gadolinium and lanthanum. Partitioning between the dissolved and particulate phases in the Rhine River and spatial propagation through the Rhine-Meuse Delta (the Netherlands), Appl. Geochem., 2014, 47, 186–197.
46. A. Slabon, S. Terweh and T. O. Hoffmann, Vertical and Lateral Variability of Suspended Sediment Transport in the Rhine River, Hydrol. Processes, 2025, 39, e70070.
47. S. N. Lane, D. R. Parsons, J. L. Best, O. Orfeo, R. A. Kostaschuk and R. J. Hardy, Causes of rapid mixing at a junction of two large rivers: Río Paraná and Río Paraguay, Argentina, J. Geophys. Res.: Earth Surf., 2008, 113, F02024.

---