A Multi-Faceted Approach to Quantifying Recovery of Stream Phytobenthos Following Acute Herbicide Incidents

: Phytobenthos, major primary producers in freshwater, may be periodically exposed to herbicides through runoff, spray drift, and leaching, but experimental work on their response to herbicides is limited. Outdoor ﬂumes were used to measure the recovery of intact phytobenthic bioﬁlms, following a short-term pulse to a glyphosate-based herbicide (GBH) and chlorotoluron (CLT), singly or as a mixture (GBHC). Two days following the removal of the herbicides, signiﬁcant declines in biomass and rates of areal gross maximum photosynthesis were recorded in GBH and GBHC-treated ﬂumes. Herbicide pulses negatively impacted the biomass of phytobenthos, measured as chlorophyll a , which in turn may have alleviated resource limitation, evidenced by signiﬁcant increases in biomass-speciﬁc rates of gross maximum photosynthesis. After 4.5 days, biomass concentrations were similar in controls and treatments, indicating bioﬁlm recovery in the GBH and GBHC bioﬁlms, though the proportion of green algae relative to diatoms was signiﬁcantly higher in GBH and GBHC-impacted bioﬁlms. Signiﬁcant increases were recorded in the Trophic Diatom Index (TDI), indicating benthic diatom development along different trajectories, following the herbicide pulse. Collectively, these results indicate greater insights into herbicide impacts on phytobenthos may be gained by working with the whole phytobenthic community.


Introduction
Implementation of the Water Framework Directive [1] has been a significant step towards ensuring the sustainability of Europe's waterbodies. By quantifying the condition of selected biota, relative to their hypothetical 'expected' condition, it is possible to identify issues responsible for degradation in waterbodies and then implement appropriate measures to restore them to at least 'good ecological status'. One of the biotic elements used for assessment as part of this process in freshwater is phytobenthos (attached algae), for which benthic diatoms [2][3][4] or the entire phytobenthic community (detailed in [5]) have been used as assessment tools, now deployed in rivers across Europe. These tools mostly focus on the major pressure of eutrophication [6]. However, rivers and streams in agricultural landscapes receive a cocktail of other stressors, including herbicides, in addition to inorganic nutrients. Herbicides may reach waters from point or diffuse sources [7], and high concentrations may result from storm runoff, leaching and spray drift into rivers and streams draining agricultural landscapes [8][9][10][11], typically subjecting benthic communities to short-term episodic pulses of toxicants or in some locations, more chronic exposures.
A few studies report on the impacts of herbicides on phytobenthos in rivers [12][13][14], but still questions remain. Outside of a controlled laboratory environment, these communities will be subject to multiple stressors, and species may vary in their responses to any one

Materials and Methods
Phytobenthos were grown in gravity-fed experimental stream channels fed directly from the Mill Stream, a tributary of the River Frome, Dorset, UK (50 • 40 48 N, 2 • 11 06 W), located on the site of the former Freshwater Biological Association River Laboratory. These channels have previously been used for experimental purposes and have been shown to support the growth of representative communities of algae and macroinvertebrates [37]. Agricultural usage in the catchment includes permanent grassland with dairy or stock rearing and cereal cultivation (Wessex Water, pers. comm). Glass slides were used to grow phytobenthic biofilms. These have been demonstrated to act as suitable substrates for colonization by microalgae [38]. Eight semi-glazed slides were located in purpose-built slide holders (phytometers), twenty of which were placed in the stream channels and left to colonise for ten days (Figure 1a).
shown to support the growth of representative communities of algae and macroinvertebrates [37]. Agricultural usage in the catchment includes permanent grassland with dairy or stock rearing and cereal cultivation (Wessex Water, pers. comm). Glass slides were used to grow phytobenthic biofilms. These have been demonstrated to act as suitable substrates for colonization by microalgae [38]. Eight semi-glazed slides were located in purpose-built slide holders (phytometers), twenty of which were placed in the stream channels and left to colonise for ten days (Figure 1a). (b) Self-contained flumes with recirculating water for short-pulse herbicide treatment.

Flume Design
It was not possible to make direct applications of herbicides into the stream channels as the river water which fed these was released directly back to the parent stream. Consequently, purpose-built self-contained flumes were set up alongside the experimental stream channels, modified from the design of Gold et al. [39]. Twenty flumes, each 1 m in length, were individually fed by a 14 L reservoir of unfiltered river water, recirculated through each channel by means of a pump (Minipond 900 m Blagdon, Dorking, UK) ( Figure 1b). Following the initial colonisation phase of phytobenthos in the stream channels, slide holders housing the colonised slides were transferred to the flumes and left to acclimate for a period of 12 h. Allocation of each slide holder and allocation of herbicide treatments between the twenty flumes were randomised. Three additional colonised slides were randomly allocated to each flume for measurement of photosynthesis and respiration rates.
After a 12 h acclimation period, herbicides were added, and biofilms were exposed for a 12 h period spanning daylight hours. Five control flumes received river water only (CON), five were dosed with glyphosate (GBH), five with chlorotoluron (CLT), and five received both herbicides in the same concentrations as the single application treatments

Flume Design
It was not possible to make direct applications of herbicides into the stream channels as the river water which fed these was released directly back to the parent stream. Consequently, purpose-built self-contained flumes were set up alongside the experimental stream channels, modified from the design of Gold et al. [39]. Twenty flumes, each 1 m in length, were individually fed by a 14 L reservoir of unfiltered river water, recirculated through each channel by means of a pump (Minipond 900 m Blagdon, Dorking, UK) ( Figure 1b). Following the initial colonisation phase of phytobenthos in the stream channels, slide holders housing the colonised slides were transferred to the flumes and left to acclimate for a period of 12 h. Allocation of each slide holder and allocation of herbicide treatments between the twenty flumes were randomised. Three additional colonised slides were randomly allocated to each flume for measurement of photosynthesis and respiration rates.
After a 12 h acclimation period, herbicides were added, and biofilms were exposed for a 12 h period spanning daylight hours. Five control flumes received river water only (CON), five were dosed with glyphosate (GBH), five with chlorotoluron (CLT), and five received both herbicides in the same concentrations as the single application treatments (GBHC). Glyphosate (Monsanto UK, Limited) was applied as the formulation Roundup ® to give a final concentration of 64.7 mg L −1 of the active ingredient (Isopropylamine salt of glyphosate), in line with the reported EC50 value determined for the green alga Raphidocelis subcapitata [40]. The formulation also contained the surfactant polyoxyethylene tallow amine (POEA). Chlorotoluron (analytical grade, PS2078, Supelco; Sigma-Aldrich Company Ltd., Gillingham, UK) was applied to give a final concentration of 8.5 µg L −1 , in line with EC50 value reported for Raphidocelis subcapitata [41]. Short-term pulses of herbicides are common when heavy rainfall follows herbicide application in river catchments [10]; hence, a 12 h exposure period represented a realistic scenario. River water from the Mill Stream was used to replace the herbicide application at the end of the exposure period, after which Phycology 2023, 3 28 further exchanges with fresh river water were made every 24 h to allow for recruitment of new cells to the phytobenthic biofilms [42].
A portable pH meter (H198185, Hanna Instruments, Woonsocket, RI, USA) was used to monitor the pH of the flume reservoir water throughout the incubation period. The average pH was 8.3 (standard error = 0.02 pH units) and did not change in any of the herbicide treatments compared to the untreated controls. Activity measurements were conducted after 36 h following the removal of the herbicides. At the end of the experiments (4.5 days after the removal of the herbicides), all phytobenthic biofilms were removed from the slides using a sterilised toothbrush. Biofilm material from eight slides within each flume was pooled, homogenised to slurries, and subsampled to determine the biomass and species composition within each flume. Figure A1 shows a schematic overview of the study design.

Biofilm Functioning (Activity Measurements)
Metabolic functioning (photosynthesis and respiration) was measured over a day, with the first measurements made 36 h post-herbicide removal. A Hansatech Chlorolab 2 system (Hansatech Instruments, Norfolk, UK) with a DW2/2 liquid-phase chamber run in dark and light conditions was used to measure rates of photosynthesis and respiration on homogenised subsamples, using samples from the additional colonised slides with a replicate taken from each set of the five flumes. A circulating water jacket around the incubation chamber maintained a constant temperature of 20 • C, the ambient temperature in the flumes. Subsamples were stirred constantly with a magnetic stirring bar. Light was supplied with LH11/2R red LEDs, and calibration of the light source was undertaken using a QRT1/PAR light sensor (Hansatech Instruments Ltd., Norfolk, UK). Field measurements of irradiance were taken beside the flumes using a Li-250 cosine-corrected probe (Li-Cor, Nebraska, USA). Ten steps of increasing irradiance were set for lab-based measurements to a maximum irradiance of 1372 µmol photons m −2 s −1 . Standard plots of gross photosynthesis versus irradiance were produced, and the theoretical coefficient of light-use efficiency (α), the maximum photosynthetic rate (P max ), and the light saturation parameter (Ik) were determined using the hyperbolic tangent model equation of Jassby and Platt [43].

Biomass Measurements
To normalise the photosynthesis and community respiration rates, sub-samples of the homogenate were used to determine algal biomass, using total chlorophyll a (Chl. a) and ash-free dry mass (AFDM).
At the end of the experiment (4.5 days), additional fresh slides were sampled for biomass measurements. For estimates of dry mass (DM), subsamples were filtered onto pre-ashed Whatman GF/F filters (47 mm) and dried at 105 • C for 24 h (n = 20, with 5 replicates per treatment). Filters were subsequently placed in a muffle furnace for 2 h at 500 • C for determination of organic content as AFDM. Chlorophyll (Chl. a) was determined by filtering a known volume of suspension through a Whatman GF/F filter. The filters were placed in a 15 mL conical tube and frozen at −20 • C to be stored and transported back to the lab for further processing. Pigments were extracted in pre-chilled acid-free 100% ethanol in an ice bath working in the dark. The absorbance was measured at 665 nm and 750 nm using a spectrophotometer (UVmini1240, Shimadzu, Kyoto, Japan). Pigment concentrations were calculated following Lorenzen [44] and adjusted for surface area of colonisation. Rates of photosynthesis were normalised to biomass using Chl. a whilst rates of community respiration were normalised using AFDM and Chl. a. Biomass measurements were used to calculate yields at the two selected times during the experimental period.

Bacterial Abundance
Bacterial abundance was determined using epifluorescence microscopy [45,46]. Subsamples were collected at the end of the experiments (4.5 days) and filtered on 0.22 µm black polycarbonate filters, stained with DAPI (Invitrogen, Life Technologies), and imaged using an Olympus BX41 epifluorescence microscope with an F-view attachment (IICCD camera (Olympus). Cells were counted using the cell image analysis software CellProfiler [47,48] and adjusted for surface area.

Algae and Cyanobacterial Species Composition
Subsamples of the pooled biofilm material were preserved in glutaraldehyde (2% final concentration). Areal-based densities were determined using a Sedgwick-Rafter Cell Counting Chamber and an inverted microscope (Wilovert, Wetzlar, Germany) following Brierley et al. [49]. Observations and measurements on uncleaned material were made at magnifications of 40×, 100×, and 400× to cover the size range of species. Where possible, identification was carried out to genus or species level using a variety of keys and guides, including Cox [50] for Bacillariophyta, whilst Gutowski and Foerster [51,52] were used for cyanobacteria and green algae. Nomenclature was adjusted to current norms using Algaebase [53]. At this level of magnification, it was impossible to identify all organisms to species, and some operational taxonomic categories were devised. Small spherical green algae (5-13 µm diameter) are described in this study as unidentified Chlorophyta, as no correct phylogenetic classification was possible using light microscopy [54][55][56]. Benthic diatoms identified for uncleaned material were assigned to one of four guilds, including motile, adnate/short-stalked, erect, and filamentous. Cell biovolume was determined according to the formulae devised by Hillebrand et al. [57]. Areal cell densities (cells cm −2 ) and biovolume (µm 3 cm −2 ) were calculated for biofilms, randomly selected from three of the five replicate flumes for each of the treatments and the controls (n = 15, with three replicates for each treatment). Based on the cell numbers and biovolume, the relative abundance and biovolume of the individual species were calculated. By grouping the taxa, total and relative cell numbers and biovolume for three different phyla, Chlorophyta, Bacillariophyta, and Cyanophyta, hereafter referred to as greens, diatoms, and cyanobacteria, were investigated.
For higher resolution, sub-samples of biofilm material were digested to determine the relative abundance of cleaned diatom valves to species where possible, and data were subsequently used to further explore differences in assemblages. Subsamples of biofilms from each of the channels (n = 20, with five replicates for each treatment) were prepared following standard methods for removal of organic matter [2]. A minimum of 300 valves were identified and their relative abundance was determined by light microscopy using a Leica DLM 500 at a magnification of 1000×. Diatoms were identified to species where possible. Some difficulty was experienced when distinguishing between some small (<10 µm in length) and often poorly silicified naviculoid diatom species using light microscopy. A category was created designated 'sNav', which included Mayamaea atomus, Sellaphora nigri (formerly Eolimna minima), Fistulifera saprophila, and other poorly resolved taxa that fell within the same size range. Difficulty in confidently assigning some of the Nitzschia to species level led to the creation of a second category, designated 'Niag', which included Nitzschia. cf. capitellata and N. cf. intermedia forms that overlapped in size and other features. The Trophic Diatom Index (TDI), an overall indication of inorganic nutrient enrichment, was determined using DARLEQ [2].

Data Analyses
Data were checked for compliance with parametric testing prior to analyses. Homogeneity of variance was tested using Levene's test, and normality of residuals was tested using the Kolmogorov-Smirnov test. Analysis of variance (ANOVA) was performed to look for differences between measures of biofilm biomass and activity. Post hoc Dunnett's tests were applied to identify differences between phytobenthos in the control flumes and treatments (p = 0.05). A non-parametric difference test (Kruskal-Wallis) was applied for data that could not be transformed. Where appropriate, Mann-Whitney (M-W) tests were carried out subsequently to identify differences between groups. Spearman's rank (SR) correlations (R) were carried out for comparisons on selected biotic variables. Tests for significance were carried out using Minitab (version 19).
The entire phytobenthos assemblage, based on cell biovolume, and a sub-set, using the relative abundance of cleaned diatom frustules as proxies for the entire assemblage [5], were ordinated using Non-metric multidimensional scaling (NMDS) with sparse species, i.e., those occurring in 4 or less samples, removed. Bray-Curtis similarity measures and MRPP were used to identify between treatment group significance. Both data sets were analysed using Pisces software Community Analysis Package (CAP), version 6.2.4, (Pisces Conservation Ltd., New Milton, Hants, UK) for NMDS and PC-ORD, version 6 (Wild Blueberry Media LLC, Corvallis, OR, United States) for MRPP, respectively. Bubble plots were used to visualise between group differences of selected taxa in the cleaned diatom sub-set using the Pisces software. For the cleaned diatom data set, TWINSPAN analysis was applied to identify species contributing to group splits and Spearman's Rank correlations were carried out between the NMDS sample scores along Axis 1 and 2 and the Trophic Diatom Index (TDI) scores [4].

Areal Biomass
Thirty-six hours after the removal of the herbicides, areal biomass measured as Chl. a, in biofilms exposed to CLT, was similar to concentrations in unamended control samples (CON) ( Table 1). Post hoc tests indicated that Chl. a concentrations were significantly lower in the biofilms treated with GBH and the combined application, GBHC, compared to controls (ANOVA, F (3,16) = 11.92, p < 0.001). Significant differences were also recorded in the dry mass (DM) of biofilms (ANOVA, F (3, 16) = 6.31, p = 0.005) 36 h into the recovery period. Biofilms from CON channels had a higher mean dry mass compared to GBH and GBHC treatments, but no difference was recorded between controls and CLT treatments (Table 1). Similarly, mean ash-free dry mass (AFDM) was significantly different between treatments (ANOVA, F (3,16) = 5.22, p = 0.012), with post hoc tests indicating that the AFDM of CON samples was significantly higher than biofilms treated with GBHC (Table 1). Table 1. Phytobenthos biomass included areal chlorophyll a (Chl. a), dry mass (DM), and ash-free dry mass (AFDM), measured 36 h after removal following a short-term (12 h) herbicide pulse compared to control samples (CON) of unamended stream water. Herbicide treatments include chlorotoluron (CLT), glyphosate (GBH), and a mixture of the two herbicides (GBHC); values are means ± SE, n = 5, * = p < 0.05.

Metabolic Responses to Herbicide Treatments
Photosynthesis versus irradiance curves for areal rates of gross photosynthesis, following a 36 h recovery period, are illustrated in Figure 2a. A significant difference was measured in the maximum rates of areal gross photosynthesis (GP max area ) derived from the plots (ANOVA, F (3,12) = 6.43, p = 0.005). Post hoc tests indicated reduced values for GP max in GBH-and GBHC-treated phytobenthos compared to the controls, with a 44% and 55% reduction in the mean GP max rates, respectively, but no differences in CLT samples relative to controls. A weak significant difference was recorded for the light utilisation coefficient (a) (ANOVA, F (3,12) = 3.52, p < 0.05), but no differences were identified between treatments and controls in post hoc Dunnett tests. No significant differences were found for the light-saturation coefficient (I k ) between controls and any treatments. plots (ANOVA, F(3,12) = 6.43, p = 0.005). Post hoc tests indicated reduced values for GPmax in GBH-and GBHC-treated phytobenthos compared to the controls, with a 44% and 55% reduction in the mean GPmax rates, respectively, but no differences in CLT samples relative to controls. A weak significant difference was recorded for the light utilisation coefficient (a) (ANOVA, F(3,12)= 3.52, p < 0.05), but no differences were identified between treatments and controls in post hoc Dunnett tests. No significant differences were found for the lightsaturation coefficient (Ik) between controls and any treatments. When normalising to biomass (Chl. a), the maximum rate of biomass-specific photosynthesis (SP max ) was about two times higher than controls with a mean ± SE of 34.4 ± 3.2 and 30.1 ± 3.7 mg O 2 mg −1 Chl. a h −1 ) in GBH and GBHC treatments, respectively (ANOVA, F (3,13) = 14.06, p < 0.001; Figure 2b).
Areal rates of community respiration were correlated to the biomass (Figure 2c; Spearman Rank correlation = 0.812, p < 0.001), but there was considerable within-group variation. Significantly higher rates were measured for CON biofilms compared to those treated with GBH and GBHC treatments when normalised to AFDM (Figure 2d; ANOVA, F (3,16) = 17.83, p < 0.001), but no differences were recorded for those treated with CLT, compared to controls. No significant differences were recorded between treatments and control when using Chl. a for normalising community respiration (p > 0.05).

Biomass of Autotrophs and Heterotrophs
After a longer recovery period, no significant differences were recorded in biomass measures, including Chl. a, dry mass, or ash-free dry mass of phytobenthos, across controls and treatments on the fourth day post-herbicide removal (p > 0.05; Table 2). The values of Chl. a and DM were lower in GBHC but still not significantly different due to within-group variation. Table 2. Areal biomass measures including chlorophyll (Chl. a), dry mass (DM), ash-free dry mass (AFDM), 4.5 days after end of exposure to a short-term (12 h) herbicide pulse compared to control treatments not exposed to herbicide; values are means ± SE, n = 5.

Enumeration of Bacterial Heterotrophs
Bacterial abundance ranged from 0.82-11.45 × 106 cells cm −2 and was variable between replicated biofilms. Overall, the total numbers of bacteria did not differ between treatments at the end of the recovery period (p > 0.05).

Taxonomic Composition of Autotrophs within Biofilms
At the end of the experiment, diatoms and green algae (belonging to the classes Trebouxiophyceae and Chlorophyceae) were the dominant autotrophs within the biofilms, with only minor representation by filamentous cyanobacteria in all treatments. In total, more than 60 different taxa were identified from the preserved material. The total biovolume for each guild and the mean relative abundance, based on biovolume, are shown for each treatment in Table A1.
No significant differences were measured between total areal cell densities (ANOVA, F (3,11) = 1.14, p = 0.39) or total areal biovolume (ANOVA, F (3,11) = 0.71, p = 0.573) across all treatments and controls for phytobenthic biofilms based on cell counts after the 4.5 day recovery period (Table A1). The ratio of the total numbers of diatoms to green algal cells was significantly lower in biofilms previously treated with GBH and with combined GBHC (ANOVA, F (3,8)  In biofilms treated with GBH and GBHC, there was a reduction in the relative biovolume of taxa assigned to the adnate/short-stalked guild and the erect guild, and an increase in the relative biovolume of motile species (Table A1). Further detailed analyses of taxonomic changes, based on cleaned diatom frustules, are reported in (Section 3.2.4).
The most common green algal taxa are listed for each treatment in Table A2. Amongst the green taxa, the chlorophyte Gongrosira sp. and the chlorophyte Desmodesmus sp., green flagellate species, and unidentified spherical chlorophytes increased in total biovolume within the GBH-and GBHC-treated biofilms compared to CON (Table A2).
A non-metric multidimensional scaling (NMDS) plot, based on the entire species composition of the phytobenthos, had a moderate level of stress (0.166). Two clusters were observed, with samples from the CON and CLT grouped close together and a second cluster comprising the GBH and GBHC assemblages located to the left of Axis 1 (Figure 4). MRPP could not confirm that differences between groups were greater than those within (A = 0.03, p = 0.227), but pairwise tests trended towards a significant difference between CON and GBHC (p = 0.11). 0.15 (SE = 0.02) in GBH and 0.13 (SE = 0.04) in combined GBHC-treated biofilms. The ratio of the total biovolume of diatoms to green algae was ~ eight times higher in the untreated controls, which had a ratio of 4.82 (SE = 2.29), and in biofilms treated with CLT, which had a ratio of 4.66 (SE = 2.06), compared to biofilms treated with GBH, with a ratio of 0.96 (SE = 0.29), or combined GBHC, with a ratio of 0.61 (SE = 0.19) (ANOVA: F(3,8) = 6.15, p < 0.05) (Figure 3b). In biofilms treated with GBH and GBHC, there was a reduction in the relative biovolume of taxa assigned to the adnate/short-stalked guild and the erect guild, and an increase in the relative biovolume of motile species (Table A2). Further detailed analyses of taxonomic changes, based on cleaned diatom frustules, are reported in (Section 3.2.4).
The most common green algal taxa are listed for each treatment in Table A3. Amongst the green taxa, the chlorophyte Gongrosira sp. and the chlorophyte Desmodesmus sp., green flagellate species, and unidentified spherical chlorophytes increased in total biovolume within the GBH-and GBHC-treated biofilms compared to CON (Table A3).
A non-metric multidimensional scaling (NMDS) plot, based on the entire species composition of the phytobenthos, had a moderate level of stress (0.166). Two clusters were observed, with samples from the CON and CLT grouped close together and a second cluster comprising the GBH and GBHC assemblages located to the left of Axis 1 (Figure 4). MRPP could not confirm that differences between groups were greater than those within (A = 0.03, p = 0.227), but pairwise tests trended towards a significant difference between CON and GBHC (p = 0.11).

Taxonomic Analysis of Pennate Diatom Fraction, Based on Cleaned Valves in Phytobenthic Biofilms
A non-metric multidimensional scaling (NMDS) plot, based on the species composition of benthic diatoms, had a low level of stress (0.111) and revealed a clear separation of the four treatment types along Axis 1. MRPP confirmed that differences between groups were greater than those within (A = 0.189, p < 0.001). Pairwise tests showed significant differences between CON and GBH (p < 0.01) and GBHC (p < 0.01), but no significant difference was identified between CON and CLT communities.
The first division of a TWINSPAN analysis (Eigenvalue (λ) = 0.162) separated four of the CON biofilm samples from all other herbicide-treated biofilms. The main indicator species leading to the split was the pooled category of small naviculoids 'sNav'. A bubble plot for this group (Figure 5a) illustrates a greater representation of this group in the biofilms treated with herbicides. The median relative abundance (%RA) of 'sNav' was significantly different between the treatment groups (Kruskal-Wallis, p = 0.002); the post hoc test showed a relatively lower median %RA of 'sNav' in CON biofilms compared to CLT-

Taxonomic Analysis of Pennate Diatom Fraction, Based on Cleaned Valves in Phytobenthic Biofilms
A non-metric multidimensional scaling (NMDS) plot, based on the species composition of benthic diatoms, had a low level of stress (0.111) and revealed a clear separation of the four treatment types along Axis 1. MRPP confirmed that differences between groups were greater than those within (A = 0.189, p < 0.001). Pairwise tests showed significant differences between CON and GBH (p < 0.01) and GBHC (p < 0.01), but no significant difference was identified between CON and CLT communities.
The first division of a TWINSPAN analysis (Eigenvalue (λ) = 0.162) separated four of the CON biofilm samples from all other herbicide-treated biofilms. The main indicator species leading to the split was the pooled category of small naviculoids 'sNav'. A bubble plot for this group (Figure 5a) illustrates a greater representation of this group in the biofilms treated with herbicides. The median relative abundance (%RA) of 'sNav' was significantly different between the treatment groups (Kruskal-Wallis, p = 0.002); the post hoc test showed a relatively lower median %RA of 'sNav' in CON biofilms compared to CLT-, GBH-, and GBHC-treated biofilms (M-W, U = 16, 15, 15, respectively; p < 0.05). There was a highly significant positive correlation (SR 0.7, p < 0.001) between NMDS Axis 1 sample scores and the %RA of the 'sNav' category. treated flumes. There was a significant difference in the pooled median %RA of three of the more well-represented adnate species, including Cocconeis lineata, Planothidium frequentissimum, and P. rostratum ('CoPl') (Kruskal-Wallis, p < 0,05). The median %RA of (~25%) in GBHC samples was significantly higher than in that of controls (~3%) (Figure  5d, M-W, U = 16, p < 0.05). The second split in the TWINSPAN dendrogram separated three of the GBHC samples from the rest of the group due to their relatively higher representation of this adnate grouping (Eigenvalue (λ) = 0.156). The upper-storey filamentous species Melosira varians was common in some samples, though distribution was patchy across all biofilms, complicating interpretation. The average %RA of the erect diatom Gomphonema parvulum was 8% in controls, but was less-well represented in herbicide treated samples.

Trophic Diatom Index in Relation to Herbicide Application
The Trophic Diatom Index (TDI) was determined for benthic diatoms in biofilms from each of the twenty flumes. The TDI scores (mean ± SD) for the GBH (63.0 ± 3.6) and GBHC (59.7 ± 6.3) treatments were significantly higher than those for the control (49.4 ± 7.6 (ANOVA, F(3,16)= 5.18, p = 0.011). A highly significant positive correlation was recorded between Axis 1 of the NMDS ordination and the scores for the Trophic Diatom Index (Pearson's r = 0.902, n = 14, p < 0.001). Nitzschia species (including N. palea, N. paleacea, and the 'Niag' group) dominated the diatom assemblage in control flumes comprising 39% of all valves. There was a significant difference in the %RA of a pooled group of N. palea and N. paleacea ('Nipp') (Kruskal-Wallis; p = 0.01), with post hoc testing recording a higher median %RA in the CON compared to the GBHC treatment (Figure 5b; M-W, U = 38, p < 0.05). No significant differences were found between biofilms in the 'Niag' group.
Achnanthidium minutissimum was also common in control samples with a median %RA of 34%, but was also relatively common in herbicide-treated biofilms. A small significant difference was recorded between treatments for the median %RA of A. minutissimum (Kruskal-Wallis, p = < 0.05), with post hoc testing indicating a higher median %RA in CON biofilms relative to GBH-treated biofilms (M-W, U = 38, p < 0.05) (Figure 5c).
In control samples, species from the genus Cocconeis (C. lineata, C. euglypta, and C. placentula) made up 9% RA collectively but were relatively more common in herbicidetreated flumes. There was a significant difference in the pooled median %RA of three of the more well-represented adnate species, including Cocconeis lineata, Planothidium frequentissimum, and P. rostratum ('CoPl') (Kruskal-Wallis, p < 0,05). The median %RA of (~25%) in GBHC samples was significantly higher than in that of controls (~3%) (Figure 5d, M-W, U = 16, p < 0.05). The second split in the TWINSPAN dendrogram separated three of the GBHC samples from the rest of the group due to their relatively higher representation of this adnate grouping (Eigenvalue (λ) = 0.156). The upper-storey filamentous species Melosira varians was common in some samples, though distribution was patchy across all biofilms, complicating interpretation. The average %RA of the erect diatom Gomphonema parvulum was 8% in controls, but was less-well represented in herbicide treated samples.

Trophic Diatom Index in Relation to Herbicide Application
The Trophic Diatom Index (TDI) was determined for benthic diatoms in biofilms from each of the twenty flumes. The TDI scores (mean ± SD) for the GBH (63.0 ± 3.6) and GBHC (59.7 ± 6.3) treatments were significantly higher than those for the control (49.4 ± 7.6 (ANOVA, F (3,16) = 5.18, p = 0.011). A highly significant positive correlation was recorded between Axis 1 of the NMDS ordination and the scores for the Trophic Diatom Index (Pearson's r = 0.902, n = 14, p < 0.001).

Impacts of Chlorotoluron on Phytobenthic Biofilms
Chlorotoluron is a phenylurea herbicide, inhibiting plastoquinone binding of the D1 protein, which plays a key role in photosystem II (PSII). It was anticipated that photosynthetic activity would be reduced in flumes treated with CLT, as the nominal concentration of CLT applied for this study (8.5 µg L −1 ) had previously been shown to bring about a 50% reduction in the growth of the chlorophyte Raphidocelis subcapitata in a single-species laboratory experiment [41]. Crucially, in standard single-species toxicity tests for algae (e.g., [58]), the species are typically exposed to the toxicant for a period of 72 or 96 h. The CLT was applied for a relatively short period (12 h) in our experiment, and this exposure period may not have been long enough to cause significant detrimental effects. However, there is evidence of a 90% uptake of herbicides by phytoplankton within the first hour of exposure to herbicides in the water column [59], giving sufficient time for the CLT to have penetrated cells exposed to the overlying water, including potential inocula, as well as more erect or larger filamentous species of diatom in the upper storey of the biofilms [26]. More recently, higher IC50-96-h values of 50 and 80 µg L −1 have been reported for the growth of a freshwater green alga Ankistrodesmus fusiformis and a marine diatom Amphora coffeaeformis [60], and Anton and Alia [61] reported an EC50-96-h of 100 µg L −1 for the freshwater green alga Chlorella pyrenoidosa. Together, these data suggest a high degree of species-specific toxicity to this toxicant. To the best of our knowledge, no directly comparable data are available for benthic freshwater species.
Whilst the thickness of the biofilms in the flumes may have acted as a barrier for penetration of toxicants such as herbicides [62], any inocula present in the phytoplankton or located in the upper regions of the biofilms (e.g., the high-profile erect and filamentous species) may have been more immediately impacted by herbicides. By contrast, cells embedded deep within an extracellular polymeric (EPS) matrix [21,63,64] may have been exposed to much lower doses. The degree of impact of toxicants may therefore be dependent on the location of an organism within the biofilm [65].
The lack of any significant differences in any of the biomass measures (Chl. a, DM, AFDM) for CLT, 4.5 days post-herbicide removal, provided further evidence that this herbicide had no measurable impact or did not impact recovery when applied to phytobenthos in a short-term pulse during this study. However, there were some shifts, albeit not significant, in the composition of the entire algal community following the addition of CLT when compared with biota in control biofilms (Figure 4). There was considerable heterogeneity in algal composition within each treatment, and no significant differences were recorded in any of the compositional measures. It is therefore recommended for further studies to increase sample sizes to increase statistical power. The total abundance of the green algae Gongrosira was reduced in biofilms treated with CLT compared to controls. This variability in phytobenthic assemblages is well recognised and is an important consideration when making assessments of ecological status in field-grown phytobenthic communities. Periodic toxicity testing on biofilms from any site over a longer time-period would better incorporate the natural spatial and temporal variability in the composition of the phytobenthic assemblages, and potentially provide greater insight regarding the 'winners' and 'losers' following herbicide exposure in specific sites.

Impacts of a GBH and GBHC on Phytobenthic Biofilms
Responses in key functional and structural metrics for GBH-and GBHC-treated biofilms were similar. Glyphosate (N-phosphonomethylglycine) affects the shikimic (or shikimate) acid pathway, impacting the enzyme involved in the production of aromatic amino acids required for protein synthesis. Negative effects of GBH formulations on plants include impacts on photosynthesis and pigment production, e.g., chlorophyll and carotenoids [66]. Our finding of a 44-55% significant reduction in areal rates of photosynthesis (GP max area ), 76-78% reduction in Chl. a, and 45-53% reduction in dry mass in GBHand GBHC-treated phytobenthic samples, respectively, within two days into the recovery period would suggest similar pathways to impact were operating in phytobenthic algae in the flumes. Roundup ® has been demonstrated to be considerably more toxic to the green alga Selenastrum capricornutum (now identified as Raphidocelis subcapitata) (EC50-96-h 5.8 mg L −1 ) and the marine diatom Skeletonema costatum (EC50-96-h 1.85 mg/L) than either the IPA salt of glyphosate or technical glyphosate [67]. In the flumes, toxic impacts were measured two days after the removal of the applied nominal concentration (64.5 mg L −1 a.i.) of Roundup ® . Such high concentrations equate more with cases of accidental spillage [68], and typically, concentrations of glyphosate in surface waters worldwide are generally much lower (~1 mg L −1 ) [28], though much higher environmental concentrations have been measured in freshwater in Argentina (105 mg L −1 ) [69]. At one extreme, Romero et al. [68] reported an EC50-96 h of 55.62 mg L −1 for the growth of a tolerant green alga Chlorella kessleri following exposure to the GBH formulation ATANOR ® . At the other end of the scale, Tajnaiová et al. [70] recently reported an IC50-72-h of 0.267 mg L −1 for Desmodesmus subspicatus using Roundup ® ClassicPro. These extreme examples serve to illustrate the wide variation in sensitivity (>200-fold difference) of these organisms, making it very challenging to predict the impacts of GBH-formulated herbicides on multi-species communities. This variation in sensitivity may not simply relate to the impacts of glyphosate per se. At least 750 different GBH formulations have been produced [71], which may contain additional (often undisclosed) adjuvants, which can induce greater toxicity than the glyphosate alone [72].
The significant increase in the biomass-specific maximum photosynthetic rates in GBHtreated samples (PS max ) two days after removing the herbicide could be accounted for by the marked reduction in biomass in the GBH-treated biofilms. Goldsborough and Brown [73] measured maximum photosynthetic rates in biofilms of relatively low biomass (2-3 µg cm −2 Chl. a), and reduced rates at relatively high biomass concentrations (~8-16 µg cm −2 ), which is comparable with our findings. The reduction in biomass following GBH addition would have alleviated stresses in the remaining under-storey biota, exposing them to an increased irradiance and potentially also removing other constraints, e.g., nutrient limitation or gas diffusion. Bacteria are known to break down glyphosate to form amino-methylphosphonic acid (AMPA) and release phosphorus [74], which could also have potentially contributed to the phosphorus pool.
During the period of herbicide exposure in flumes, it is possible that organisms could have left the biofilms in response to unfavourable conditions. A drift response may enable organisms to escape from resource-limited conditions [75] or unfavourable habitats [76]. The biofilms lacked a cohesive physiognomy, increasing susceptibility to sloughing events [77]. Sloughing of phytobenthic material has been observed following the application of titanium dioxide nanoparticles to phytobenthos in experimental flumes grown from the Mill Stream, Dorset, UK [78]. Evidence suggests that under-storey algae in phytobenthic biofilms may retain photosynthetic potential despite being potentially shaded as the biofilms mature, and reduction in over-storey biofilm may result in an improved photosynthetic performance [79].
Significantly higher areal community respiration rates were measured in control biofilms compared to GBH and GBHC herbicide-treated phytobenthos circa two days into the recovery period. Respiration rates in GBH and GBHC treatments were only significantly lower when normalised using organic content rather than Chl. a. Whilst bacterial numbers were not significantly different between treatments at the end of the experiment, it is possible that exposure to these herbicide treatments may have, at least in the period during and immediately after removal, affected the functioning (e.g., secondary production) of the bacteria community. We did not examine the heterotrophic protist component of the biofilms, but these may also have been negatively impacted. However, a much greater sensitivity to the Roundup ® formulation has previously been demonstrated in algae, in comparison to bacteria and protozoa [67].

Changes in Species Composition of Phytobenthos
Differences in the diatom species composition, visualised by NMDS, 4.5 days after removal of the herbicides provided evidence of species-specific responses, with greater spatial separation for GBH-and GBHC-treated phytobenthos from controls than for CLT. The significantly higher proportion of three adnate diatom species in GBHC-treated biofilms indicated a marked difference in the diatom assemblage compared to other treatments after the 4.5-day recovery period (Figure 5d), suggesting exposure of under-storey species in these biofilms and potentially increased growth and/or new recruitment of these adnate species. Given that fresh river water was provided for all flumes, new inocula could have rapidly colonised the biofilms, or more tolerant species may have proliferated, replacing herbicide-sensitive species. Diatoms of the genus Gomphonema were common in control biofilms but were sparse or absent from most GBH-and GBHC-treated biofilms. Wood et al. [80] tested the impact of glyphosate on a homogenised riverine benthic diatom assemblage using indoor laboratory microcosms. Based on the percentage of healthy cells, they recorded an EC50-48 h for glyphosate of 5.3 µg L −1 and >500 µg L −1 for Gomphonema spp. (a mixture of G. parvulum and G. minutum) but an EC50-48 h of >500 µg L −1 (the highest concentration applied) was recorded for G. gracile, demonstrating different within-species tolerance. Navicula cryptotonella, the only naviculoid present in their river community, was the most tolerant genus, corroborating our observations that other species in this genus were tolerant to GBH and GBHC.
Green algae colonised the biofilms rapidly as spaces in the biofilms became available, as shown by the marked decline in the ratio of diatoms to green algal cells in GBH and GBHC biofilms. Other researchers have also recorded changes in the composition and/or biomass of phytoplankton and phytobenthos following herbicide applications including GBH [81][82][83]. Here, we demonstrate evidence of a full recovery, at least for biomass measures, with biomass concentrations matching those of controls after the 4.5-day recovery period, which implied rapid growth of remaining taxa or rapid re-colonisation from the inoculum over a very short space of time in the GBH and GBHC biofilms. Our results concur with those reported by Debenest et al. [12], who stated that whilst no overall effects of biomass may be measured in the long-term following herbicide application, it is possible that increases in the biomass of certain taxa (e.g., green algae) may compensate for inhibition in others (e.g., diatoms). Here, some unicellular Chlorella-like green algae and Desmodesmus species opportunists, recognised as taxa capable of rapidly colonising newly exposed areas of biofilms [22], increased in abundance. The timeframe would have allowed for several divisions of tolerant cells already in the biofilm. The recirculating nature of the flumes used for our experiments, and the daily replacement of stream water, allowed also for the possibility of both emigration and re-immigration of biota in addition to new colonists, as would occur under natural conditions in streams after a peak event of herbicide runoff, affording greater environmental realism.
Although much research has been carried out on the impacts of glyphosate and GBH on phytoplankton and phytobenthos communities, there is still no clear message regarding the impacts of glyphosate-containing herbicides. Negative, positive, or no impacts have been reported, with outcomes likely to depend on the composition and state of the cohort at the start of the experiment (e.g., the health of cells and development stage of the community, the composition of the GBH, nutrient composition of the freshwater, length of exposure to the toxicant, and experimental design (indoor, outdoor, flowing, or static) replenishment of new inocula [27,72,81,82,[84][85][86][87][88]. Incorporating an experimental setup where fresh river water is replaced during the recovery period, as utilised here, adds greater environmental realism for measuring the potential for recovery of herbicide-impacted biofilms.

Implications for Ecological Assessment
Whilst we show that phytobenthic biofilms had the capacity to make a rapid recovery from an acute toxicity event, we also provide evidence (through TDI measurements) that changes observed in species composition may link to an enhanced growth of species that are tolerant to higher nutrient concentrations, given the breakdown product of glyphosate could lead to the release of phosphorus [74,89]. The fact that no measurable growth was recorded in CON-and CLT-amended biofilms over the 'recovery' period compared to biofilms treated with GBH or GBHC alluded to a possible limitation to growth in the former and alleviation of this constraint in the latter. The reduction in the relative abundance of the nutrient-sensitive diatom Achnanthidium minutissimum in GBH-treated flumes and the significant increase in the relative abundance of nutrient-tolerant small naviculoids resulted in a significant increase in the TDI in GBH-and GBHC-treated biofilms, although we cannot rule out that these species may just be more tolerant to glyphosate, or less sensitive to any associated adjuvants. In the 2-year period prior to this experiment, improvements in phosphorus stripping led to a reduction in the mean soluble reactive phosphorus (SRP) concentrations (49 µg L −1 ) in the river Frome, falling below the target of 60 µg L −1 set by the Environment Agency (EA) for chalk rivers [90,91]; though, the possibility of nutrient limitation within thicker biofilms in our experiments could not totally be ruled out.
Based on observations made on river biofilms situated in catchments of intensive agricultural activity, Roubeix et al. [92] commented that it is difficult to separate the impacts of nutrients from those of herbicides as the two are highly correlated. We concur with this view, as it was not possible to confirm if the increased biomass in the GBH and GBHC treatment in the recovery period was a direct result of (i) herbicide action per se possibly favouring increased growth of herbicide-tolerant species with reduced competition from more sensitive species; (ii) phosphorus release from glyphosate; (iii) greater availability of nutrients already in the water column due to reduced biofilm thickness; or (iv) alleviation of other growth constraints generated from the reduction in biofilm thickness (e.g., light) [93]. Relative to control biofilms, the ratio of the total biovolume of diatoms to green algae in GBH and GBHC biofilms was about 8-fold lower, showing a significant shift away from a diatom-dominated biofilm (Figure 3). It is possible that green algae were more tolerant to herbicides, and they may just fill the gaps left by emigrating diatoms or sloughed biofilm. Measurements of significant changes in biomass and metabolic activity in the short-term contrasted with the lack of difference measured in biomass a few days after application of the herbicides, demonstrating the rapid recolonization within impacted biofilms. The evidence presented here suggests that the timing of measurements when assessing impacts to herbicides may be critical for determining impacts. Furthermore, the selection of parameters to measure change is critical. Had we only used biomass measures, such as Chl. a, DM, and AFDM, measured 4.5 days into the recovery period, the marked changes in the species composition would have been missed.
The challenge for ecological assessment is to simultaneously evaluate the combined impact of these multiple stressors and to isolate the roles of individual stressors in order to give catchment managers realistic guidance on appropriate remediation measures and enable the development of specific protection goal options, as also formulated by the tasked European Food Safety Authority [94]. This study complements our previous work (e.g., [2]) by demonstrating how the impact of herbicides on benthic communities may confound the interpretation of diatom indices designed to reflect the predominant nutrient-organic gradient in European streams. It also raises questions about the feasibility of disentangling herbicide impacts from those of other pressures using the current suite of methods, which focus largely on taxonomic composition [3].
Not only do herbicides cause different responses dependent upon their modes of delivery and action, but we may also anticipate a broader range of outcomes in reality compared to what has been reported from single-species tests or diatom-only approaches. However, even when using our more integrative design looking at green algae, diatoms, and cyanobacteria, factors that we held as constants (pre-exposure time, duration of exposure, and herbicide concentrations) will all be highly variable in the field. We studied two herbicides, whereas river communities will be exposed to many agrochemicals, often in mixtures, and the pathways which transfer these from the site of application to streams will also influence the biota in other ways (rainfall will alter nutrient concentrations and increase hydrological stress on biofilms, and potentially increase sediment loads). Whilst our results suggest that herbicide impacts will be reflected by the current UK assessment method [2], we are pessimistic about the potential of taxonomic-based methods generally to provide a strong diagnostic capability. Our observation of a strong shift in proportional dominance between diatoms and green algae in response to herbicides suggests a need for a broader taxonomic range in ecological assessment than the current narrow focus on diatoms [95]. Our results on photosynthetic activity suggest that functional aspects should be an important part of future toxicity testing and indicate, contrary to what might have been expected, that productivity may be enhanced in the short term by the application of certain herbicides, and changes in growth may have important consequences for higher trophic guilds [96]. The provision of river-grown biota when testing herbicide impacts offers more realistic conditions for quantifying impacts at a community level, compared to laboratory experiments with a single species.
An alternative view of the situation is that all spatial studies of algae and cyanobacteria in regions of intensive agriculture will almost certainly include intermittent effects of herbicides as one of several factors contributing to the principal stressor gradients. Such gradients may be correlated with nutrients. Herbicides are rarely included as variables in such studies due to the difficulties in chemical measurement. However, their influence will be manifest in the biota, as we demonstrate here. The ability of the biota to detect pressures missed by conventional chemical analyses is often quoted as a reason in favour of ecological assessment, a point that underlies the rationale for the WFD specifying using phytobenthos. The utility of ecological assessment methods is often demonstrated by the strength of correlation or association with chemical gradients [97]. Ironically, herbicides may well confound the relationship between the biota and variables such as phosphorus which are widely used to calibrate or validate methods. The long-term solution may be to recognise the limitations of community-based approaches as diagnostic tools and, instead, to focus on their capacity to provide rapid overviews of the condition of the aquatic biota [3]. Funding: We are very grateful to the Freshwater Biological Association (FBA) for providing 50% of the funds to support a PhD Studentship (Grant #SM6505). Grateful thanks are also given to the Bristol Centre for Agricultural Innovation (BCAI) (formerly Long Ashton Research Station, LESARS) for providing 50% of the funds towards a PhD Studentship (Grant #HF1743).

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Acknowledgments:
We would like to thank Lutz Kupferschläger for designing the graphical abstract and Figure A1. We are very grateful to the Freshwater Biological Association, (FBA) for the use of their site facilities and the support of the River Lab staff at the FBA River Lab, Dorset.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of this study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

Appendix A
Phycology 2023, 2, FOR PEER REVIEW 18 Appendix A Figure A1. Schematic illustration of the study design. Figure A1. Schematic illustration of the study design.