This study investigated the Holtemme, a 47 km long river originating in the Harz mountains of Saxony-Anhalt, Germany. The river is considered near-pristine in the mountainous region, while downstream anthropogenic impact increases due to wastewater discharge, rectification, settlements, and widespread farming activities (Inostroza et al., 2016; Karthe et al., 2017; Wollschläger et al., 2017). In the mountainous region, low concentrations of caffeine and the insecticide diethyltoluamide indicate a mild impact by hikers (Weitere et al., 2021). The sediment in the Holtemme is mostly sandy with pebbles and cobbles (Kunz et al., 2017). Five small streams flow into the Holtemme in the study area, all of which contribute little to the water volume of the river. During the sampling period in the summer season of 2022, a stream entering the Holtemme between Site 5 and 6 had dried up completely. For this study, two sampling sites were selected in the “pristine” part of the river, one directly at the discharge point of an activated sludge-based, tertiary WWTP, one about 150 m downstream of the WWTP, and two others 8 and 13 km further downstream. The latter site is at the upstream outskirts of the next city (Halberstadt), which has its own WWTP. All six sampling sites have been previously described in detail (Haenelt et al., 2023).

Planktonic and biofilm samples were collected during a dry period in the summer season of 2022 on five sampling days (August 9th, three consecutive days from August 23rd to 25th, and September 1st) at all six sites. Rain events were rare and the average rainfall was less than 2.3 mm per day (Supplementary Figure 1). Water temperature and pH were measured at all sites using a SenTix41 probe (Xylem Analytics, Germany). NO₃-N, total phosphorous, and water depth were obtained from the MOBICOS stations located at Sites 1 and 5 (Fink et al., 2020), and physicochemical parameters from Site 3 (WWTP discharge) were kindly provided by the WWTP Silstedt.¹

Planktonic samples were collected in sterilized and pre-rinsed 1-L glass flasks at an approximate water depth of 5 cm. Biofilm samples were collected with a spatula from the water-facing side of stones covering the riverbed and transferred into tubes containing 1 mL 1xPBS. To obtain a sufficient amount of biomass, approximately 10 cm² of surface area was used for biofilm sampling at Sites 1, 2, 5, and 6. At Sites 3 and 4, which were characterized by dense biofilms, approximately 1 cm² was sampled. All samples were stored in a thermobox with cold packs and transported to the laboratory within 3 h. Immediately upon arrival, biofilm samples were stored at 8°C until DNA extraction within the next 2 days. Planktonic samples were further treated as follows. One liter of river water from each site except Site 3 (500 mL) was filtered using 0.22 μm pore size PSE PALL filters (diameter 47 mm, Pall Corporation, New York, NY, USA) and stored at −20°C until DNA extraction. Additionally, 100 mL flow-through was collected for each site and stored at 8°C overnight before solid phase extraction (SPE) and sulfonamide concentration measurements via HPLC-MS.

Flow-through was concentrated 100-fold using SPE with Oasis HLB 6cc 500 mg columns (Waters, CT USA) before HPLC-MS was carried out as described in Haenelt et al. (2023). In brief, a calibration curve from 0.01 μg/L to 2 μg/L of a sulfonamide mixture containing sulfamethoxazole (SMX), sulfadiazine (SDZ), and sulfamethazine (SMZ) was generated, and all samples were measured in technical triplicates. HPLC-MS analysis was carried out using a Zorbax Eclipse Plus Rapid Resolution HT-C18 (100 mm × 3 mm, 1.8 μm) column on a 1,260 Infinity II HPLC (Agilent Technologies, Santa Clara, CA, USA) coupled to a QTRAO 6500 + MS/MS (AB Sciex, UK).

The biofilm samples were vortexed and centrifuged for 5 min at 13,000 g to harvest microorganisms. The supernatant was discarded and the sediment was resuspended in 100 μL BE buffer (Macherey Nagel, Germany). For planktonic samples, filters were cut into small pieces and divided into two parallel 2-mL tubes containing 200 μL BE buffer and approximately 0.15 g of 1-mm zirconium beads (Carl Roth, Germany) each. The individual tubes were agitated for 20 s at 4 m/s on a FastPrep 24 instrument (MP Biomedicals, Germany) and centrifuged in a benchtop centrifuge (neoLab, Germany). The supernatant from parallel tubes was collected and pooled. To maximize the final DNA yield, an additional 100 μL of BE buffer was added to the tubes containing the filter pieces, vortexed, and centrifuged again, and the supernatant was pooled with the previous one. DNA extraction was performed using the NucleoSpin Microbial DNA Kit (Macherey Nagel, Germany) following the manufacturer’s protocol.

The diversity of bacterial communities in planktonic and biofilm samples was determined for all sites and sampling campaigns by preparing a 16S rRNA gene amplicon library as described (Haenelt et al., 2023). We conducted a nested PCR of the variable region V3 using the primer pairs specified in Supplementary Table 1. Performing a nested PCR was necessary to obtain sufficient DNA yield for the planktonic samples. Amplicons were sequenced with the NextSeq 500/550 High Output Kit v2.5 on an Illumina NextSeq 550 instrument. Bacterial community analysis was performed using QIIME 2 (Bolyen et al., 2019), DADA2 (Callahan et al., 2016), and the taxonomic classifier “Silva 138 99% OTU full-length sequences” (Bokulich et al., 2018; Robeson et al., 2021). Phyloseq (McMurdie and Holmes, 2013) was used for further data analysis in R. The dataset for Site 1 on August 9th was removed due to low sequencing depth. Rarefaction curves can be found in Supplementary Figure 2. For statistical analyses of alpha diversity measurements, we performed Kruskal-Wallis and Dunn’s tests using the package rstatix (Kassambara, 2023). An overview of the statistical significance of differences in alpha diversity, both between sampling sites and within each sampling site, comparing planktonic and biofilm samples, can be found in Supplementary Table 2. To compare microbial beta diversity in biofilm and planktonic samples, we applied non-metric multidimensional scaling (NMDS) of Bray-Curtis dissimilarities using the vegan package (Oksanen et al., 2022). To assess seasonal influences on bacterial community structure and diversity in the surface water compartment, we did NMDS with data obtained for planktonic samples in this study and data from our previous study from the 2020/2021 winter season (Haenelt et al., 2023). PERMANOVA tests were performed to evaluate the statistical significance of differences in beta diversity between sampling sites, within each sampling site comparing both surface water and biofilm compartments (Supplementary Table 3), as well as between the summer and winter season.

Absolute abundances of sul1, sul2, intI1, and the 16S rRNA gene were determined with SYBR Green-based quantitative real-time PCR using established primers (Supplementary Table 1) and a protocol from Haenelt et al. (2023). For planktonic samples, copy numbers (CN) were normalized to CN/100 mL river water. For biofilm samples, CN were normalized to CN/sample and are therefore not directly comparable with each other. Relative abundances of sul1, sul2, and intI1 were calculated by dividing their absolute abundances by those of the 16S rRNA gene. The statistical significance of differences in absolute and relative ARG abundances between sampling sites and seasons was determined using the package rstatix (Kassambara, 2023). We performed statistical analyses employing Kruskal-Wallis and Dunn’s tests using the package rstatix (Kassambara, 2023). Results for the relative abundance of the target genes, encompassing comparisons of both sampling sites and riverine compartments within each sampling site can be found in Supplementary Table 4. Additionally, we compared the abundance of all four target genes quantified in the current study with data acquired during a winter season in our previous study (Haenelt et al., 2023). We aimed to identify notable differences in both absolute and relative abundance between the two seasons for each sampling site using Dunn’s test, and the outcome of these analyses is depicted in Supplementary Table 5.

Our monitoring study showed that WWTP effluent was the primary source of SMX, nitrogen, and phosphorus pollution. Downstream of the discharge point, our data show decreasing SMX as well as nutrient concentrations. The highest SMX concentration was found at Site 3, with a median concentration of 20 ng/L, resulting in 17.7 ng/L SMX at Site 4. This finding is in line with previous studies on the Holtemme (Krauss et al., 2019; Beckers et al., 2020; Švara et al., 2021) and with numerous publications which showed that WWTP effluent is the primary source of antibiotic pollution in aquatic ecosystems (Carvalho and Santos, 2016; Rodriguez-Mozaz et al., 2020; Zhuang et al., 2021; Larsson and Flach, 2022). Sulfonamides are typically found at concentrations below 100 ng/L in anthropogenically impacted river waters (Li et al., 2022). In our previous study during the 2020/2021 winter season (Haenelt et al., 2023), the median SMX concentration at Site 3 was 8.6 ng/L and therefore lower than during the summer season. Lower antibiotic concentrations in the effluent during the winter season could be explained by dilution effects due to precipitation in the WWTP catchment (Hong et al., 2018). Taking the SMX concentration as an indicator of the proportion of wastewater in the river, WWTP effluent would account for nearly 90% of the river water at Site 4 in the summer season of 2022. This value is similar to the value of 75% estimated on the basis of the flow rates provided by the State Office for Flood Protection and Water Management Saxony-Anhalt.² Downstream of the discharge point, SMX concentration decreased. However, the decrease was not statistically significant and the decline rate of SMX concentration (11% from Site 3 to Site 6) was lower than during the winter season (73%), presumably because of high water evaporation together with the lack of precipitation. SDZ and SMZ concentrations did not exceed 1.3 ng/L at any of the sampling sites (Figure 1A). SMX and SMZ concentrations at Site 4 are consistent with a previous study conducted on the Holtemme (Tousova et al., 2017). There, samples were taken approximately 1 km downstream from our Site 4, the mean SMX concentration was 8 ng/L and SMZ concentrations were below the limit of quantification. SDZ concentrations were previously not quantified.

During the hot and dry summer season in 2022, water depth was on average below 10 cm in the Holtemme. Water temperature upstream of the WWTP varied between 13 and 16°C, while wastewater discharge at Site 3 (21°C) resulted in an elevated water temperature of 19°C at Site 4. Downstream at Site 5 (17°C), water temperature increased further, peaking at 22°C at Site 6. The pH in the river changed from pH 7.3 at Site 1 to pH 8.8 at Site 6. Nitrate concentrations were 1.2 mg/L upstream of the WWTP, 4 mg/L in the wastewater effluent, and 2.8 mg/L downstream at Site 5. These values are in the lower range of previously reported nitrate concentrations along the river stretch (Brase et al., 2018). Total phosphorous showed a similar pattern with low concentrations of 0.02 mg/L at Site 1 and an increase to 0.1 mg/L at Site 5 after WWTP discharge of 1 mg/L (Figure 1B).

The bacterial community structures of planktonic and biofilm samples are reported individually, using Shannon index and species richness based on the number of observed Amplicon Sequence Variants (ASVs) to assess bacterial diversity (Figure 2). In planktonic samples, both indices showed the same trend. Highest median values of 6.4 and 1,379 were found at Site 1, whereas lowest values were found at Site 3 (median 5.5 and 1,113). The decrease in Shannon index from Site 1 and 2 to Site 3 was statistically significant (Dunn’s test: p ≤ 0.01 for Site 1, p ≤ 0.05 for Site 2). Downstream of the WWTP, both indices increased steadily, reaching values of 6 and 1,354 at Site 6. The steady increase of both indices downstream the discharge point indicates a recovery of the riverine bacterial community, which was demonstrated in our previous winter monitoring study on the Holtemme using a microbiome recovery model (Haenelt et al., 2023). In contrast to planktonic samples, both indices decreased in the biofilm samples from Site 1 (5.8 and 1,038) to Site 2 (5 and 782). Values at Site 3 (5.7 and 1,157) and nearby Site 4 (5.9 and 1,334) are similar to those at near-pristine Site 1, while further downstream the WWTP, both indices again decreased in a statistically significant manner (Dunn’s test: p ≤ 0.05) with medians of 5.3 and 803 at Site 6. Notably, the biofilm communities at Sites 1, 2, 5, and 6 showed significantly lower alpha diversity in comparison to the corresponding planktonic samples (Dunn’s test: p ≤ 0.05). The number of observed ASVs at Site 1 was the only exception, showing a non-significant difference. Conversely, polluted Sites 3 and 4 show much more uniform alpha diversity in both planktonic and biofilm-forming bacterial communities. We hypothesize that the increase in both indices in the surface water compartment downstream from the WWTP might be attributed to biases during rarefaction of the sequencing data. The high abundance of certain taxa in wastewater may have prevented detection of rare taxa, which were found further downstream in the river with higher evenness (Cameron et al., 2021).

The composition of bacterial communities along the river stretch was impacted by anthropogenic pollution. This is evident in the dissimilarities of beta diversity across different sites, as shown in the NMDS plot in Figure 2C. Influencing factors could be foremost the discharge of wastewater, flow rate, physicochemical parameters, and rectification (Luo et al., 2020). A clear trend could be observed with Sites 1 and 2 on one side and WWTP effluent and nearby Site 4 on the other side of the plot, whereas Sites 5 and 6 showed very similar community structure and are positioned in between the upstream and WWTP-impacted sites. Significantly distinct beta diversity was observed in both planktonic and biofilm-forming microbial communities across all sampling sites, except when comparing Site 3 to Site 4 (PERMANOVA: p ≤ 0.05). NMDS analysis also revealed a statistically significant separation between planktonic and biofilm-forming bacterial communities at all sampling sites (PERMANOVA: p ≤ 0.01). However, the change in community composition seemed to be driven more by the anthropogenic pollution rather than the type of compartment, i.e., surface water vs. biofilm.

The taxonomic compositions of the samples corresponded to the degree of anthropogenic influence along the river. The 50 most abundant phylotypes for each site and their corresponding relative abundances can be found in Supplementary Table 6. Biofilm and planktonic communities at Site 1 were dominated by members of the classes Alphaproteobacteria, Gammaproteobacteria, Clostridia, and Bacilli. In the biofilm samples, approximately 4% of the high-quality reads were derived from Cyanobacteria, whereas in the planktonic samples, Cyanobacteria were not among the 50 most abundant phylotypes (>0.26% of total read abundance). At the family level, bacteria typically exhibiting an aerobic, heterotrophic metabolism predominated in both compartments (e.g., Chitinophagaceae, Comamonadaceae, and Xanthobacteraceae). Noteworthy, we found an ASV affiliated with Escherichia-Shigella in every sample from the surface water compartment at Site 1 (average read abundance of 0.3%). Since identical or similar ASVs were not found in any of the other 55 samples from our campaign, it seems unlikely that they are sampling/laboratory contaminants in the planktonic samples from Site 1. While it remains to be elucidated to which genus the Escherichia-Shigella ASV belonged, and whether it was derived from an autochthonous population, the members of that genus could be a source of the occasionally detected sul1, sul2, and intI1 genes at the site (see below). At Site 2 just downstream the city of Wernigerode, the communities of the surface water and biofilm compartments included several phylotypes often associated with warm-blooded animals (e.g., various Ruminococcaceae) in addition to those found at Site 1. The WWTP effluent at Site 3 was dominated by phylotypes of the family Peptostreptococcaceae and the order Lactobacillales in planktonic samples, while in biofilm samples most abundant phylotypes were associated with Peptostreptococcaceae, Nitrospiraceae, and Rhodobacteraceae. Downstream the discharge point, Peptostreptococcaceae remained a dominant family, but Sphingomonadaceae and Xenococcaceae were also detected with increasing abundance. Arcobacteraceae (phylum Campilobacterota) became highly abundant in biofilm samples (Figure 3). Members of this family have been identified as predominate hosts of sul1, sul2, intI1, and probably class 1 integrons in sewage-impacted water (Hultman et al., 2018; Knecht et al., 2023). Overall, the relative abundance of families with known hosts of integrons slightly decreased from 17% in both biofilm and surface water compartment at Site 1 to an average relative abundance of 12% downstream the WWTP, thus potential hosts for ARGs were already present upstream the WWTP in an area with little anthropogenic impact. Some of these families such as Aeromonadaceae and Rhodocyclaceae have been implicated as prevalent carriers of the indicator genes in and downstream of wastewater treatment systems (Zhang et al., 2018; Lo et al., 2022; Knecht et al., 2023).

Differences in beta diversity of planktonic communities from this sampling campaign and those analyzed in our study of the winter season 2020/2021 (Haenelt et al., 2023) are shown in Figure 4. Bacterial communities in summer and winter season were distinct, with statistically significant differences at all sampling sites (PERMANOVA: p ≤ 0.05). Summer communities cluster more than those acquired in winter season, which is most likely related to the shorter sampling period. This becomes evident in the spatial arrangement of data points in the NMDS plot depending on the sampling date in the winter season, where samples were collected over approximately 3 months. The summer sampling campaign lasted only about 1 month. In contrast to upstream sites, the bacterial community in the WWTP effluent were more stable over time and in between seasons, as data points are located close to each other in the plot. As WWTP performances are often stable over time, one would expect the bacterial community structure to be distributed more closely around an average, as previously shown by Wang et al. (2021). At downstream Sites 4, 5, and 6, site-specific distances are much smaller than the time-specific differences. These findings are consistent with a monitoring study conducted on a river in southern Germany (Reichert et al., 2021), illustrating the substantial impact of seasonality on bacterial community structure and diversity.

The anthropogenic pollution pattern in the Holtemme is evident from the prevalence of sul1, sul2, intI1, and the 16S rRNA gene in planktonic samples. All four genes showed a similar pattern with low CN or CN below limit of quantification at Sites 1 and 2, highest CN in wastewater and decreasing CN downstream of the WWTP (Supplementary Figure 3). The increase in absolute abundance was statistically significant for sul1, sul2 and 16S rRNA gene from upstream Sites 1 and 2 to heavily wastewater-impacted Sites 3 and 4 (Dunn’s test: p ≤ 0.01). The median CN of 16S rRNA gene at Site 6 was only 16% of the CN at Site 4. Similar declines were observed for sul1, sul2, and intI1 (12, 10, and 6%, respectively). The absolute CN of the 16S rRNA gene was higher at the wastewater discharge point (median 4.1 × 10⁷ CN/100 mL) than in the river water with a median of 4.3 × 10⁶ CN/100 mL at Site 5. In contrast to our previous monitoring study in the winter season 2020/2021, the absolute abundance of all four target genes in summer season declined more rapidly from Site 4 to 5. This observation is of special interest, as we would have expected the opposite, a faster decrease in the absolute abundance during winter season, due to dilution effects. Absolute abundance for biofilm samples is not depicted because a normalization to CN/sample does not allow for direct comparison between different sampling sites. In the summer season, relative abundance against 16S rRNA gene showed a similar pattern for all three target genes sul1, sul2, and intI1 in planktonic samples, with the highest relative abundance at Site 3 (median 2.8, 2.6, and 1.3%, respectively) and decreasing relative abundance from Site 4 to 5 (Figure 5). Wastewater discharge at Site 3 resulted in significantly increased relative abundances of sul1 and sul2 when compared to both upstream sites (Dunn’s test: p ≤ 0.05). The decrease in absolute and relative abundance downstream the WWTP is most likely attributed to dilution effects, but additional removal mechanisms like sedimentation, cell death by predation or inactivation by UV light might also occur (Lee et al., 2021).

Although both, the relative and absolute abundances of these genes are decreasing, the median values remain largely unchanged 8 km (Site 5) and 13 km (Site 6) downstream the WWTP. These observations show that ARGs can be persistent in the aquatic environment even far downstream the point source of pollution (>10 km). Sabri et al. (2020) measured the abundance of several ARGs in a Dutch river and absolute abundances of sul1, sul2, and intI1 gene copy numbers remained constant over 20 km of the river stretch regardless of the season. These findings differ from our previous monitoring study in winter season 2020/2021, where absolute abundances of sul1 and sul2 decreased to copy numbers below the limit of quantification at Site 6, 13 km downstream the WWTP. Potential causes for differing observations could be the fraction of wastewater discharge, and the amount of precipitation and additional inflow from tributaries depending on the season. In addition, the absolute abundance of the 16S rRNA and sulfonamide resistance genes was significantly higher in the WWTP effluent and at all sampling sites downstream the WWTP in the summer compared to the winter season (Dunn’s test: p ≤ 0.05). This increase is most likely a consequence of higher abundance of microorganisms in the WWTP during summer season. The relative abundance of sul1 and sul2 was also higher in the summer season at Site 3 and downstream the WWTP. For sul2 the difference was statistically significant (Dunn’s test: p ≤ 0.05). We estimate an increase in the proportion of WWTP effluent at Site 4 from ∼27% in winter to about 75% in the summer season (see above). Therefore, seasonal differences like the increased proportion of wastewater and thus ARG-carrying wastewater-derived microorganisms in the river increased the abundance of sulfonamide resistance genes. In the face of climate change with less average precipitation in summer seasons in Europe, one must assume that the proportion of wastewater in rivers will increase even further (European Environment Agency, 2021), leading to an increase of ARGs in rivers.

In biofilm samples, similar relative abundances were observed for sul1 and sul2, however, at higher relative abundance at Site 3 (both genes at 3%) and Site 4 (2 and 4.7%) and a more rapid drop from Site 4 to Site 5 (both 0.3%). The relative abundance of both genes was significantly increased at Site 3 and 4 when compared to the pristine Site 1 (Dunn’s test: p ≤ 0.01). The elevated relative abundance in vicinity to the discharge point indicates that ARG-carrying microorganisms were mainly accumulated in biofilm directly at the WWTP, whereas those remaining in surface waters over longer distance did not or only barely enter biofilms at a later stage. Similar findings were reported by Guo et al. (2018), where the authors found higher ARG abundances and detection frequencies in biofilm compared to planktonic samples. The accumulation of antibiotics in the biofilm compartment, in combination with dense microbial growth and accumulation of mobile genetic elements might render biofilms an important reservoir for ARGs downstream WWTPs (Gillings et al., 2009). However, the spatially restricted occurrence of sulfonamide resistance genes sul1 and sul2 in biofilms in close vicinity to the discharge point suggests that this was not the case in our study.

In contrast to the sulfonamide resistance genes, intI1 was detectable upstream and downstream of the WWTP without significant differences in the relative abundance between sampling sites. However, wastewater discharge increased its relative abundance from a median of 0.9% at Site 2 to 1.5% at Site 4. All anthropogenically impacted ecosystems have become polluted with the clinical class 1 integron in recent years (Gillings, 2018). For example, Cacace et al. (2019) conducted a study on the abundance of nine ARGs and the intI1 gene in 16 WWTP effluents from 10 European countries. Similar to the biofilm samples we studied, they found a background level of intI1 resistance already present upstream of the WWTP, which further increased when wastewater entered the river. Possible reasons for the occurrence upstream of a WWTP could be anthropogenic impact from the city of Wernigerode or leisure activities at the river (Weitere et al., 2021). The prevalence of the intI1 gene in biofilm throughout the river stretch is of concern if correlating with the abundance of class 1 integron, as the integron could potentially carry various ARGs as gene cassettes, and biofilms provide suitable conditions for horizontal gene transfer (Gillings et al., 2008; Chaturvedi et al., 2021).

During our monitoring study, we were limited in information regarding flow rates for all sampling sites, which would have been needed to calculate the loads of ARGs in CN/day as previously done by Lee et al. (2021). This would have further enabled us to determine the environmental fate of ARGs (e.g., dilution of wastewater in the river vs. additional removal mechanism) more precisely. In addition, single-cell approaches such as epicPCR or directGeneFISH (Moraru et al., 2010; Spencer et al., 2016) could be implemented in future studies to link the host phylotypes with the ARGs type and to better assess the risk of ARG transfer to human pathogens.