The response of microphytobenthos to physical disturbance, herbicide, and titanium dioxide nanoparticle exposure

The microphytobenthos that form transient biofilms are important primary producers in intertidal, depositional habitats, yet we have only a limited understanding of how they respond to the cumulative impacts of the growing range of anthropogenic stressors to which they are exposed. We know even less about how the temporal alignment of exposure – such as duration and exposure sequence – may affect the response. Estuarine biofilms were cultured in mesocosms and exposed to the herbicide glyphosate and titanium dioxide (TiO 2 ) nanoparticles in different sequences (glyphosate-first or TiO 2 -first), as well as in the presence and absence of physical disturbance. We found that at environmentally realistic chemical concentrations


Introduction
Transient biofilms, formed by an assemblage of photosynthetic unicellular eukaryotes and prokaryotes, and termed microphytobenthos (MPB), are a critical ecological component of unvegetated intertidal flats (Hope et al., 2019;Pinckney, 2018).A significant proportion of ecosystem primary production for many estuarine systems is undertaken by the MPB (MacIntyre and Cullen, 1995), which are often dominated by epipelic diatoms (Hubas et al., 2018;Paterson and Hagerthey, 2001) and cyanobacteria (Watermann et al., 1999).The MPB are often the basis of the estuarine food chain and an important food source for deposit feeders such as Corophium volutator and Hediste diversicolor (Gerdol and Hughes, 1994;Smith et al., 1996), which themselves serve as prey to predators such as shorebirds (Peer et al., 1986).The MPB also play a key role in the biological stabilisation of cohesive marine sediments through the exudation of extracellular polymeric substances (EPS) (Chen et al., 2017b;Paterson, 1989;Tolhurst et al., 2002).
Estuarine systems are prone to climate change associated variables such as sea level rise and weather extremes (Fujii, 2012;Wetz and Yoskowitz, 2013), and these pressures are compounded by other anthropogenic stressors such as the discharge of pollutants, plastics and xenobiotic compounds (la Farré et al., 2008;Macleod et al., 2016;Solomon and Thompson, 2003;Tovar-Sánchez et al., 2013), physical disturbance such as bait digging and boating (Beukema, 1995;Walker et al., 1989;Wynberg and Branch, 1997), and eutrophication linked to excess nutrient loading (Anderson et al., 2002;Howarth et al., 2011).Any such stimulation that causes stress to a population or ecosystem can be regarded as a stressor, and these can be further categorised as natural or anthropogenic.The type and complexity of stressors found in nature has increased dramatically since the industrial revolutionexamples are climate change, novel pollutants such as pharmaceutical and personal care products (PPCPs), and anthropogenic noise (Solan and Whitely, 2016).Whilst the relevance of different factors will vary from system to system, fundamentally there are likely to be multiple drivers of change operating simultaneously.However, our understanding of stressor impacts on marine ecological systems has until recent years been dominated by single stressor studies (Crain et al., 2008;Gunderson et al., 2016), which is problematic in terms of predicting how systems will respond to stress, as non-additive synergistic and antagonistic interactions between multiple co-occurring stressors are difficult to anticipate.Increasingly this knowledge gap is being filled, although more work on intertidal systems is needed (Alestra and Schiel, 2015;Kenworthy et al., 2016;Vye et al., 2017).
Whilst the number of multiple stressor experiments involving intertidal biofilms are limited, freshwater diatom-containing biofilms have been the subject of studies which have examined the effects of stressors such as nutrient enrichment, sediment loading, herbicides and light intensity on community structure and diversity (Magbanua et al., 2013;Wagenhoff et al., 2013;Winkworth et al., 2015).However, studies on the impacts of stressors on freshwater habitats cannot necessarily be extrapolated to marine and intertidal systems.Not only do communities vary between these systems, but freshwater and marine biofilms appear to have differing functional attributes (Spears et al., 2008).For example, sediment stabilisation in marine systems takes place by increasing the critical shear stress for erosion foremost at the sediment surface (Tolhurst et al., 2003;Underwood and Paterson, 1993) while the stabilisation of freshwater sediments becoming increasingly more significant with depth below the flocculent surface layers (Gerbersdorf et al., 2008(Gerbersdorf et al., , 2009)).
In addition to the lack of studies on the cumulative effects of multiple stressors on intertidal biofilms, few studies incorporate more than a pair of stressors simultaneously despite synergistic interactions being more common in experiments with more than two stressors (Crain et al., 2008), or consider temporal factors such as the likely scenario of asynchronous stressor exposure (Gunderson et al., 2016); this is despite the potential for these factors to influence cumulative stressor impacts (Breitburg et al., 1998;Brooks and Crowe, 2019;Fukami, 2001).This may in part be a result of the difficulty in carrying out fully factorial experiments with more than two stressors, due to the rapid increase in treatments required to test each interaction combination.The number of treatments increases further if additional conditions, such as the exposure sequence of stressors or intensity fluctuations, are incorporated into the experimental design, not to mention the need for a sufficient level of replication to maintain sufficient statistical power.This complexity does not mean that such approaches should be ignored but does serve to make such efforts more difficult.
Mesocosm studies, whilst necessarily limited in their degree of environmental realism, can be a useful tool to gain a basic understanding of the magnitude of stressor impacts, to identify potential interactions, and test hypotheses in a controlled setting (Hancock et al., 2021;Vallino, 2000).In this study, intertidal MPB biofilms were grown on sediment in a mesocosm system and exposed to the herbicide glyphosate and titanium dioxide nanoparticles (TiO 2 ) in different sequences, to test for an effect of ordering on cumulative stressor impacts.A tertiary stressor, physical disturbance, was also introduced and applied to half of the treatments, to determine whether an additional source of stress would influence the cumulative effects of the first two stressors.Increased sedimentation has been shown to interact cumulatively with glyphosate to have varying effects on different freshwater algal guilds (Magbanua et al., 2013), and resuspension of glyphosate, which bonds strongly with sediment particles (Vereecken, 2005), may cause it to become bioavailable again after an initial exposure pulse.Glyphosate was selected as a stressor due to its long-term, widespread use as a commercial herbicide (Duke and Powles, 2008) with moderate marine persistence (Mercurio et al., 2014), to which microalgae (including benthic diatoms) are sensitive (Magbanua et al., 2013;Tsui and Chu, 2003;Wood et al., 2016a).Titanium dioxide was also tested as a representative metal-oxide that leeches into the aquatic environment in the form of nanoparticles, through its use in commercial pigments (Kaegi et al., 2008;Sharma, 2009) and cosmetics such as sunscreens (Tovar-Sánchez et al., 2013).There is evidence that microalgae, including diatoms, are negatively impacted by exposure to metal-oxide nanoparticles (Baker et al., 2014;Deng et al., 2017).Disturbance as a stressor was selected to imitate the physical disturbance caused by human activities such as cockle raking (Cowie et al., 2000;Kenworthy et al., 2016) or trampling (Wynberg and Branch, 1997).Stressor effects were quantified by determining a number of functional attributes related to MPB: sediment stability; photosynthetic capacity; and chlorophyll content.
In aqueous environments, glyphosate and TiO 2 nanoparticles have been shown to interact through adsorption of the herbicide onto the metal nanoparticles, and the enhanced degradation of glyphosate to aminomethylphosphonic acid (AMPA) (Ilina et al., 2017).AMPA itself may be of a similar environmental toxicity to glyphosate, but degrades eventually to inorganic phosphate, ammonium, and carbon dioxide (Battaglin et al., 2014), potentially rendering it less toxic to the MPB over time.The exposure sequence may therefore affect glyphosate bioavailability, for example, if it is adsorbed by TiO 2 on application rather than having the direct toxic action expected and reducing the detrimental effect to the MPB.The diatom species assemblage was also quantified at the end of the experiment, as these communities have been shown to vary in response to potential stressors such as pH, light and nutrients (Lange et al., 2011;MacDougall et al., 2017).The objective of this work was to investigate how the microphytobenthos assemblage might be impacted by different exposure scenarios, and to consider whether these scenarios warrant further investigation under more environmentally realistic conditions.Because TiO 2 NP exposure was expected to have a negative effect on MPB, but combination with glyphosate was expected to reduce the potency of the herbicide, the interaction would be identified by the presence or absence of a significant sequencing effect under each of the physical disturbance scenariosif TiO 2 rendered glyphosate less potent more rapidly than under natural degradation alone, net impacts would be weaker in the treatments in which TiO 2 was first added to the mesocosms.The specific hypotheses for these experiments were: 1. Stressor exposure would have a negative impact on biofilm condition and sediment stability as compared with a control group.2. Biofilm photosynthetic functioning and diatom species composition would vary dependant on the sequence of exposure to the stressors.3. The response to exposure sequence would be contingent on whether physical disturbance took place (an interaction between disturbance and sequence).

Mesocosm establishment
Sediment cores (20 cm diameter by 10 cm depth) were taken from the Tay estuary (Fife, Scotland) at Tayport Heath (56 • 26 ′ 26.4"N, 2 • 51 ′ 46.5"W) during low tide on 22/10/2018.Cores were wet sieved through a 0.5 mm mesh to exclude macrofauna, and the remaining sediment homogenised before being added to 10.5 cm × 10.5 cm square mesocosms that were filled to a depth of 5 cm.After a short settling period, a piece of bubble wrap (10.5 cm × 10.5 cm) was carefully placed on the sediment surface, in order to minimise disturbance and resuspension of particles upon the addition of 605 cm 3 of seawater with a salinity of 25 (equivalent to the salinity recorded at the Tayport Heath field site).The seawater was treated by passing through a sub-sand filter and a settling tank, and combined with distilled water to reach the appropriate salinity.Mesocosms were then placed in a temperaturecontrolled environment at 13 • C in a rectangular array and exposed to a 12-h artificial light cycle (06:00-18:00, PAR ~ 80 μmol m − 2 s − 1 at sediment surface) to match the approximate average seasonal conditions at the field site.
Each mesocosm was equipped with an aeration tube to gently oxygenate the water; bubbles were channelled (5 cm horizontally along the mesocosm wall at 2 cm depth) to induce some circular flow in the overlying water.After a 7-day establishment period, mesocosms were drained by siphoning and pipetting, with care taken to minimise the disturbance of the surface, in order to remove excess nutrients associated with sediment sieving and handling (Hale et al., 2014).The drained water was replaced with 605 cm 3 of fresh, filtered seawater with a salinity of 25.This water was not replaced during the experiment to avoid removing chemical stressors or any derivatives (such as AMPA) from the mesocosm after application.The mesocosms were then allowed to settle for 48 h before exposure to experimental treatments, to allow time for the epipelic diatoms to migrate to the sediment surface or settle from suspension (Consalvey et al., 2004).The lack of a tidal cycle and a constant photoperiod is likely to have attenuated the natural vertical migration of the diatom species within the mesocosms, however migratory species were still observed at the end of the experiment in each treatment and the control (e.g.species of the genera Nitzschia and Navicula, see Supplementary Materials 1) and vertical migration is somewhat persistent even in the absence of a tidal regime (Palmer and Round, 1967).The presence of a permanent depth of water was also expected to have a limiting effect on rates of MPB photosynthesis, particularly compounded by the limited PAR (photosynthetically active radiation) emitted by artificial lighting in comparison with natural illumination.As a result, stressor impacts were treated as effects relative to controls, rather than predicted effect sizes under the conditions of a more natural system, indicative of the likely direction of responses and the size of a disturbance required to induce an effect.Fifteen mesocosms were established, such that each treatment and the control were replicated N = 3 times each, with an effective N = 6 for each exposure sequence and disturbance scenario due to replication across these categoriesi.e. 3 replications of each sequence in both disturbance categories, and 3 replications of each disturbance category in both sequences (Fig. 1).

Stressor application
Selected stressors were TiO 2 nanoparticles at a concentration of 0.5 μg L − 1 in the overlying water, glyphosate at a concentration of 200 μg L − 1 in the overlying water, with or without physical disturbance.The TiO 2 concentration was based on realistic environmental concentrations (Nowack and Mueller, 2008;Sharma, 2009) as Titanium (IV) oxide anatase nanopowder (batch average 21 nm primary particle size by TEM, >99.5 % trace metals basis) supplied by Sigma-Aldrich.Glyphosate was obtained as Gallup® Biograde 360 Glyphosate weed killer and was diluted with 1 mL of seawater of 25 salinity giving a concentration of 200 μg L − 1 , a value selected to be just above the measured water residue concentration found in a tributary 2 h after overspraying with glyphosate (Feng et al., 1990).These environmentally realistic concentrations were selected to ensure there would be scope for stressors to act in concert without inducing a fully lethal effect.For glyphosate, 200 μg L − 1 has been shown to have a slightly negative effect on diatom cell health (Wood et al., 2016b), and for TiO 2 0.5 μg L − 1 was significantly lower than concentrations of 1.0-5.0mg L − 1 which have been shown to have toxic effects on the marine diatom Thalassiosira pseudonana (Galletti et al., 2016).In the mesocosms, disturbance was caused by raking channels in the sediment surface to a depth of approximately 10-20 mm, across both the horizontal and vertical plane of the mesocosm surface (~30 channels raked in each plane).Four treatments and a control were established (Fig. 1) with combinations of these stressors to test for both the effect of ordering on cumulative effects, and for the effect of adding a third stressor (physical disturbance).
On day 9, after 7 d establishment and 48 h settling, the first chemical stressor was applied to each treatment (Fig. 1), and the second applied after a further 48 h.Differences among treatments immediately after stressor application may be attributable in part to this temporal gap, but this is negligible by the end of experiment (25 days after application) as initial glyphosate effects manifest rapidlywithin 24 h (De María et al., 2006) and can be expected to impact unicellular algae more rapidly than vascular plants, certainly acting within 48 h (Wood et al., 2016a).
Glyphosate was pipetted directly into each mesocosm as a 1 mL dose, whereas titanium dioxide was added from a stock suspension in deionised water before pipetting to achieve the appropriate concentrations in the volume of water in each mesocosm.During and following pipetting, the overlying water was gently swirled with the pipette tip to evenly disperse the viscous titanium solution.Then, 48 h after the last stressor application, physical disturbance was applied to mesocosms A1 and A2 (Fig. 1).

Photosynthesis metrics
A Monitoring Pen MP 100-E (Photon Systems Instruments) was used to carry out pulse amplitude modulation (PAM) fluorescence measurements to measure biofilm photosynthetic performance at three time points after stressor exposure -3, 18 and 25 days.The overlying water within each mesocosm was removed by siphon and pipette, and photosynthetic metrics were obtained through the generation of rapid light curves of the relative electron transport rate following a 15-min dark adaption period; (see Baker, 2008;Maxwell and Johnson, 2000;Pierre et al., 1998;Ralph and Gademann, 2005;White et  were α (the photosynthetic efficiency in the light-limited region of the rapid light curve, i.e., the gradient of the light curve, relative units), E k (the light saturation coefficient μmol m − 2 s − 1 ), F v /F m (maximum quantum efficiency of PSII of dark adapted samples, unitless), and rETR max (maximum relative electron transport rate, relative units).See Consalvey et al. (2005) for more information on the derivation of these metrics.The seawater removed was subsequently returned to each respective mesocosm.

Diatom community structure
On day 25, a disposable pipette was used to collect 5 subsamples from the top 2 mm of sediment of each mesocosm, which were then combined and preserved in a 2.5 % glutaraldehyde solution in seawater.The preserved sediment was resuspended by vortexing, and 0.3 mL of the supernatant was removed into a pyrex test-tube containing 0.7 mL of distilled water, and 1 mL of saturated potassium permanganate solution was then added to break down organic matter.The solution was left for 24 h after which 1 mL of 32 % hydrochloric was added, and the solutions were slowly warmed in a laboratory oven to 70 • C and left until the solution became a clear straw-colour.Solutions were then centrifuged for 10 min at 1500 rpm, and the supernatant replaced with distilled water.This process was repeated 5 times to remove any remaining acid.The cleaned pellet was then resuspended in distilled water to form a slightly cloudy solution (approx 3 mL) which was then sonicated and passed through a 90 μm mesh to remove larger sediment particles.Thirty μl of the solution was pipetted onto a microscope cover slip and left to dry overnight.Permanent slides for microscopy were then mounted in 'Naphrax'.Three hundred diatom valves were counted per slide using the Utermöhl method (Utermöhl, 1958).Species identification accuracy was verified by reference to identification keys, type specimens held at the Natural History Museum of London and scanning electron microscope images taken at the Sediment Ecology Research Group (SERG) laboratory.

Sediment stability and chlorophyll content
The tip and bottom 0.5 cm of a 20 mL syringe (2 cm diameter) was cut away, sanded, and washed in distilled water to serve as a corer to preserve the limited surface area within the mesocosm.The corer was inserted approximately 2 cm into the sediment, before being gently rotated and lifted to remove a sediment core which was then frozen with liquid nitrogen.The upper 5 mm of the frozen cores were sectioned off using a circular diamond saw, and then weighed and freeze-dried.Chlorophyll extraction was performed following the HIMOM protocols for pigment extraction (Brockmann et al., 2004) with the following modifications: quartz sand was not mixed with the sediment-solvent mixture, and an ultrasound bath was not used.Instead, the sample was vortexed at 3 time pointsimmediately following acetone addition, after 24 h, and after 48 h.Chlorophyll a content was measured using a spectrophotometer (Brockmann et al., 2004) and expressed in terms of μg g − 1 of dry sediment.
After sediment cores were extracted, sediment surface stability was measured using a cohesive strength meter (CSM) (Paterson, 1989).The CSM determines sediment erodibility by firing a jet of seawater from within a water-filled chamber at the sediment surface and measuring the attenuation of an infrared light beam caused by sediment resuspension within the chamber.The critical eroding pressure (kPa) is considered to have been reached when the transmission of the infrared beam drops by 10 % and is expressed as the critical suspension threshold (N m − 2 ) using the equation given in Tolhurst et al. (1999).The most sensitive CSM programme "Fine 1" was used, and these measurements were performed after all others due to the physical disruption of the sediment caused by measurements.

Statistics
Statistical analysis was carried out in R (R Development Core Team, 2018).The effects of stressor sequencing and physical disturbance on biofilm condition were tested using Bayesian generalised linear mixedeffects models (Stan statistical modeling software accessed through R with the RStan and rstanarm packages) (Goodrich et al., 2018;Stan Development Team, 2018).Additional details on the Bayesian model fitting approach, selection of priors, and tests for interactions are given in Supplementary materials 1. Models were constructed with PAMmetrics as response variables, and sampling time point, physical disturbance and sequencing as fixed effects.Mesocosm replicate identities were set as random intercept effects to account for multiple observations of each mesocosm through time.Sequencing and disturbance were categorised as two-level factor parameters (TiO 2 -first or second and physical disturbance -presence or absence, respectively), and leaveone-out cross validation was used to test for the presence of interactions between disturbance and exposure sequence.
Chlorophyll a and sediment stability were modelled as general linear models without a temporal predictor, as there were no repeated measures over time, with sequence and disturbance as predictors for each variable.All other model checks were carried out as described in Supplementary materials 1. Directional effect probabilities (Pd) were calculated as the percentage of posterior samples in which the estimated effect had the same sign (positive or negative) as the median effect.Treatment effects compared against the control baseline were analysed using Tukey's Honest Significant Difference (HSD) method, where the test was applied to analysis of variance models of the response against each treatment, and at each time point for repeated measures.Diatom community composition equivalence at the end of the experiment as expressed by Bray-Curtis dissimilarity was tested by permutation multivariate analysis of variance (PERMANOVA) and visualised using non-metric multidimensional scaling (nMDS) ordination.Both methods were implemented in R using the Vegan package (Oksanen et al., 2019).

Effects of disturbance and stressor exposure sequence
The effects of disturbance and sequence, and change over time, on each functional attribute is given alongside the probability that the effect on each attribute lies in the same direction (positive or negative) as the median (Pd) (Table 1).The effect of sequence could be represented with either glyphosate or TiO 2 -first (Pd is unaffected), but the median effect size is given for TiO 2 -first relative to glyphosate-first to more intuitively show whether the presence of the nanoparticles added before the herbicide resulted in a reduced herbicide effect.Each attribute, and the control condition, is discussed in more detail below.Cross-validation did not support the presence of an interaction between physical disturbance and the stressor exposure sequence on any of the measured responses, therefore an interaction column is not included (Table 1).

Chlorophyll
Mean sediment chlorophyll a content differed more between the control and treatments than among the treatments themselves, although showed similar levels of variation across mesocosm replicates (Fig. 2).A Tukey HSD test showed the control to have a significantly greater chlorophyll content than treatments A1 (p = 0.034), B1 (p = 0.037), and B2 (p = 0.026) by almost 4 of μg g − 1 dry sediment, but not significantly greater than A2 (p = 0.169).
Among the treatments there was an 83.4 % probability that exposure to physical disturbance was associated with decreased sediment chlorophyll content, with a most probable increase of 0.73 μg chlorophyll per gram of sediment, and a 67.2 % probability that TiO 2 -first was associated with a lower chlorophyll content than glyphosate-first, with a most probable decrease of 0.35 μg chlorophyll per gram of sediment (Table 1).

Critical erosion threshold
The mean critical erosion threshold was greater in each treatment than the control, however variation in the threshold among replicates was often quite high, and a Tukey HSD test showed that there was no strong evidence for a difference between the control and any individual treatment (p > 0.1 in each pairwise comparison) (Fig. 3).
Among the treatments there was a 55.0 % probability that physical disturbance had a positive effect on the critical erosion threshold, with a most probable associated increase of 0.004 N m − 2 , and an 80.0 % probability that the TiO 2 -first sequence was associated with a greater critical erosion threshold, with a most probable associated increase of 0.027 N m − 2 (Table 1).

Photosynthetic rate in the light-limited region of the rapid light curve (α)
α tended to increase over time and was greatest 25 days after initial stressor application in each treatment except the control (Fig. 4).This was reversed in the control, with α decreasing over time.Variation in the response also tended to decrease over time, tending to be greatest 3 days after stressor application in treatments A1, A2, and B1.Tukey HSD tests provided no strong evidence for differences among the treatments or between the treatments or control (p > 0.1 in each pairwise comparison) within any given time interval.
Among the treatments there was a 72.5 % probability that physical disturbance had a negative effect on α, with a most probable decrease of − 0.006 relative units, and a 69.0 % probability that the TiO 2 -first sequence was associated with a positive effect on α, with a most probable increase of 0.005 relative units.There was a 99.8 % probability that α increased over time, with a most probable increase of 0.020 relative units per time interval (Table 1).

Maximum relative electron transport rate (rETR max )
The trajectory of rETR max over time following stressor exposure varied among treatments and the control, tending to decrease in A1, B1 and B2, but sharply increasing between 3-and 25-days post stressor exposure in A2 and the control (Fig. 5).Response variance across

Table 1
Estimated effects of physical disturbance and exposure sequence on MPB response measurements when controlling for the other respective effect.Time was also controlled and estimated as an effect between periods for repeated measurements.The intercept median is the most probable measurement value in the absence of physical disturbance, and with TiO 2 as the second stressor.Pd is the probability that the effect direction was the same as the median.mesocosm replicates also tended to increase as mean rETR max increased.Tukey HSD tests provided some evidence of a difference in rETR max between the control and treatments A1 (p < 0.01) and B1 (p < 0.05) 3 days after stressor exposure, otherwise there were no clear differences between or among treatments at any time interval (p > 0.1 s in each pairwise comparison).Among the treatments there was an 80.0 % probability that physical disturbance was associated with an increase in rETR max when controlling for exposure sequence, with a most probable increase of 2.12 relative units, and a 90.5 % probability that the TiO 2 first sequence was associated with an increase in rETR max , with a most probable increase of 3.25 relative units.There was an 86.6 % probability that rETR max decreased over time, with a most probable decrease of − 1.63 relative units (Table 1).

Light saturation coefficient (E k )
E k tended to decrease over time among the treatments, being greatest 3 days after stressor exposure, however this was reversed in the control, with E k increasing over time (Fig. 6).Tukey HSD tests provided some evidence of a difference between E k in the control and treatment A1 3 days after stressor exposure (p < 0.01), otherwise there were no clear differences between or among treatments at any time interval (p > 0.1 in each pairwise comparison, except A2 vs control at 3 days p = 0.08).
Among the treatments there was an 89.8 % probability that physical disturbance was associated with an increase in E k , with a most probable increase of 17.3 μmol m − 2 s − 1 , and a 78.7 % probability that the TiO 2first exposure sequence was associated with an increase in E k , with a most probable increase of 10.9 μmol m − 2 s − 1 .There was a >99.9 % probability that E k decreased over time, with a most probable decrease of 30.1 μmol m − 2 s − 1 per time interval (Table 1).

Maximum quantum efficiency of PSII (F v /F m )
F v /F m tended to increase over time, being greatest 25 days after stressor application, however there was variation in this parameter among replicate mesocosm at 3 days in treatments A1, A2, and B1.F v / F m was similar in the control at each time point and showed relatively limited variation among replicates (Fig. 7), however Tukey HSD tests conducted at each time interval provided little evidence for differences among treatments or between treatments and the control (p > 0.1 for each pairwise comparison and time point, except control vs A1, B1 and B2 at 25 days where p = 0.10, 0.06, and 0.09 respectively).
Among the treatments there was a 77.2 % probability that physical disturbance was associated with a decrease in F v /F m , with a most probable decrease of − 0.02, and a 64 % probability that the TiO 2 -first sequence was associated with an increase in F v /F m , with a most probable increase of 0.01.There was a 99.9 % probability that F v /F m increased over time after stressor exposure, with a most probable increase of 0.05 (Table 1).

Community composition
Across the treatments and control, 81 different diatom species were identified (Supplementary materials 2).There was variation in diatom community structure among the treatments, with the greatest variation occurring between treatments and the control (Fig. 4).Following verification of homogeneity of multivariate dispersion by ANOVA (p = 0.84), community structure was tested for centroid and dispersion equivalence by PERMANOVA and found to vary significantly (pseudo-F (4,14) = 1.57, p < 0.01, permutations = 999).Therefore, stressor application appears to have driven a significant change in the community structure (Fig. 4).

Discussion
The results of these experiments showed that non-specific stressor exposure scenarios had several different impacts on the cultivated biofilms as compared to the control.Notably, the chlorophyll a content of sediment in all treatments after 25 days was lower than the control, indicating that the MPB biomass was reduced due to stressor exposure, or that growth rate was reduced relative to the controls mesocosms.The reduced chlorophyll a content of treatments was significantly lower than   the controls in all treatments except A2, with a reductive effect of approximately 30-40 % less chlorophyll a per gram of sediment in the stressed treatments.However, there was no time point at which the maximum quantum efficiency of PSII, F v /F m , was found to differ among any treatment.This metric is routinely used as an indicator of plant health (Maxwell and Johnson, 2000) and may cautiously be used for this purpose in assessing the health of benthic diatom biofilms (Consalvey et al., 2005), thus given the notable difference in sediment chlorophyll a content between the treatments and control it might be expected that F v /F m would be greatest for the control biofilms.An explanation could be that the stressors drove a change in the biofilm community composition towards more tolerant species.However, whilst there was some evidence of a greater difference in community structure between the control and treatments than among the treatments, this was limited (no significant pairwise differences).Moreover, the five most common diatom species averaged across replicate mesocosms were identical and constituted 40 % or more of the species present, differing somewhat only in their relative frequencies.It is unclear whether this limited community shift would significantly affect measurements of fluorescence or significantly alter the biological traits of the community.However, despite detecting no pairwise significant differences by Tukey's HSD test, mean F v /F m was lower in all treatments 3 days after stressor addition than after 25 (Fig. 7) and showed greater variation among replicates, particular in treatments A1, A2 and B1, indicating a recovery/change following exposure, albeit inconsistently.Across the treatments, there was stronger evidence to suggest that F v /F m increased over time, with a 99.9 % probability that the metric was greater 25 days after exposure than after 3 when controlling for sequence and physical disturbance.
Similarly, there was a 99.8 % probability that α (Table 1) increased over time in the treatment groups but remained more stable in the controls, indicating that within days of stressor exposure there may have been negative impacts on the capacity for light harvesting under the experimentally-limited light conditions and consequent reduced rates of photosynthesis as a result of changes to cell architecture and/or pigmentation (Platt and Jassby, 1976).Though the impact may be weak, the recovery of F v /F m and α in stressed treatments relative to the controls is consistent with potential photosystem damage.The maximum relative rate of electron transport 3 days after stressor exposure was consistently higher than in the control group, particularly in treatments A1 and B1 (Fig. 6), and only dropped slightly over time.Being derived from these metrics, the light saturation parameter E k showed a similar relationship with the treatments as rETR max , and the inverse of α, although statistical evidence for a difference at any time point between the control and treatments was borderline (Fig. 5).A reduction in E k relative to the controls could indicate an impact on the capacity of the light harvesting complex.Where E k and rETR max were significantly greater in the stressed treatments at 3 days relative to the control, this suggests that soon after exposure that these biofilms still had the potential for optimal photosynthetic activity, but this would eventually decline, suggesting a delayedbut limitednegative response to the applied stressors.
To summarise: • After 25 days, sediment chlorophyll a content was lower following stressor exposure • Stressor exposure induced a limited shift in the diatom community structure • There was some evidence of damage to diatom photosynthetic apparatus, but not an immediate reduction in the maximum relative electron transport rate or E k Given that there was a compelling difference in the sediment concentration of chlorophyll a of stressed versus non-stressed treatments but only limited changes in diatom community structure and fluorescence indicators of photosystem damage, either or both factors may explain the reduction in chlorophyll.However the stressor impacts manifested, there was no strong evidence of a major consequent effect on the biostabilisation of the sediments (Chen et al., 2017a;Paterson et al., 2018) by the MPB communitiesthough the mean critical erosion threshold of the control was lower than for each treatment (Fig. 3).It is interesting to note that the stressed treatments consistently displayed greater critical erosion thresholds which would be consistent with the protective exudation of EPS (Flemming and Wingender, 2010;Steele et al., 2014), warranting further investigation as critical erosion threshold measurements tend to be noisy and require greater replication than other responses.
Unlike the effect of stressor application against the control, the chlorophyll a content of sediment at the end of the experiment varied little with stressor exposure sequence; there was only a 67.2 % probability (Table 1) that TiO 2 -first was associated with lower chlorophyll a than glyphosate-first (50 % represents complete uncertainty in, or no effect of, exposure sequence).Moreover, of the PAM fluorometry metrics only the maximum relative electron transport rate had a >90 % probability of differing between the exposure sequences, estimated to be 3 relative units greater under the TiO 2 -first sequence when accounting for physical disturbance and repeated measures.F v /F m , E k and α were also estimated to be slightly greater (or less negative) when exposed to titanium dioxide before glyphosate (Table 1), but these effects were too weak to confidently conclude that the initial nanoparticle exposure reduced the effect of the subsequent herbicide application on the MPB.Similarly, whilst there was some evidence that treatments exposed to TiO 2 -first were more resistant to erosion, the probability of the effect being positive was still well below certainty (=80 %).
The estimated effect of physical disturbance, in addition to exposure, to the chemical stressors was of a similar magnitude, and similar uncertainty, to the estimated effect of varying the exposure sequence.The strongest effects of physical disturbance were on chlorophyll a, rETR max , and E k , which were greater when exposed to physical disturbance.It is not clear, therefore, that physical disturbance was acting negatively as a stressor under these experimental conditions.This is perhaps not surprising, as the impacts of physical disturbance on a benthic system may be complex.The interface between sediment and overlying water is not discrete and despite the disruption of the biofilm, MPB in nature may be resuspended whilst maintaining productivity (MacIntyre and Cullen, 1996).In the context of a closed mesocosm any suspended algae will remain within the system and resettle.Moreover, physical disturbance of sediments may act akin to bioturbation and result in the release of nutrients (Biles et al., 2002;Caliman et al., 2011).This emphasises that the definition of a stressor is not rigid, but a function of the nature of the stressor exposure and the environmental context.There was also no strong evidence for an interaction between the effect of disturbance and the sequence of stressor exposure for any of the measurements or derived fluorescence metrics, suggesting that any effect of altering the order of chemical stressor exposure was not contingent on whether physical disturbance had taken place.As both the effects of disturbance and changing exposure sequence were weak, there is limited scope for an interaction to manifest.

Conclusions
The environmentally realistic stressor scenarios used in this series of mesocosm experiments have been shown together to elicit changes in intertidal biofilms, with consequent detrimental impacts on biomass as proxied by sediment chlorophyll a content.As the key primary producers of the unvegetated intertidal flats (MacIntyre et al., 1996;Miller et al., 1996), this has potential ramifications for grazers such as Corophium volutator and Hediste diversicolor, and their predators, as well as nutrients cycling and biologically-mediated sediment stabilisation (Hope et al., 2019).However, it is not clear that the level of physical disturbance induced in this experiment was detrimental and even had some positive effects.Experiments which altered the sequence of exposure to the herbicide glyphosate and TiO 2 nanoparticles did not demonstrate with sufficient confidence that exposure to the nanoparticles first mitigated against the effects of subsequent glyphosate exposure.Further, fully factorial, experiments with these two stressors may shed more light on possible interactions between them, and other cumulative effect manifestations at greater stressor intensities.It is also possible that the effect of a stressor is only felt if it occurs at a key moment in the life cycle of the benthic organism in question.In the case of diatoms, this may correspond to the time of sexual reproduction or to the seasons during which the supply of nutrients and/or light is more important and certain functions/genes are reactivated.Indeed, cell division has been linked to vertical migration as the optimal conditions for mitosis are both depth and light dependent (Saburova and Polikarpov, 2003), suggesting a strong potential for cell cycle synchronisation.In this case, it is quite possible that the effect of the sequence could be significant if the disturbances do not occur on a linear time scale (as is the case in this study) but rather extemporaneously in relation to a given need or function.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 1 .
Fig.1.Experimental analysis structure, visualising the stressor sequencing and physical disturbance treatment codes (A1, A2, B1, B2).Treatments A1 and A2 were exposed to physical disturbance, and treatments A1 and B1 were exposed to TiO 2 first, whereas A2 and B2 were exposed to glyphosate first.Replicates: N = 3 per treatment, effective N = 6 per exposure sequence and disturbance category.

Fig. 2 .
Fig. 2. Mean chlorophyll a content (μg g − 1 dry sediment) of sediment by treatment and days after stressor exposure with standard deviation bars.

Fig. 3 .Fig. 4 .
Fig. 3. Mean critical erosion threshold (N m − 2 ) of sediment by treatment and days after stressor exposure with standard deviations bars.

Fig. 5 .
Fig. 5. Mean maximum relative electron transport rate by treatment and days after stressor exposure with standard deviation bars.

Fig. 6 .
Fig. 6.Mean E k by treatment and days after stressor exposure with standard deviation bars.

Fig. 7 .
Fig. 7. Mean F v /F m by treatment and days after stressor exposure with standard deviation bars.

Fig. 8 .
Fig. 8. NMDS ordination of diatom community structure grouped by treatment (k = 2, stress = 0.16), where treatment codes are defined in Fig. 1.Coloured polygons connect each replicate (n = 3), such that overlapping polygons indicate similar community structures and larger polygons indicate greater dispersion within the respective treatment community.