Deciphering the origin of riverine phytoplankton using in situ chlorophyll sensors

Riverine algal groups with distinct life histories can generate unique patterns of structural and functional behavior. As such, novel methods to discriminate between these groups can improve the understanding of river ecosystem processes. We examined benthic vs. planktonic contributions to suspended algal biomass by monitoring suspended chlorophyll concentration and turbidity during 48 storm events at 2 locations with contrasting hydraulic storage associated with low‐head dams. Upstream from the dams, chlorophyll hysteresis showed concentrating effects and counterclockwise rotation, suggesting stormflow concentrated algae from benthic sources. When autotrophic conditions (P/R > 1) preceded storms, chlorophyll hysteresis switched to more proximal benthic sources (faster mobilization). Downstream of the dams, hysteresis showed greater dilution effects and more proximal sources of planktonic algae than at the upstream site, in contrast to high similarly in turbidity hysteresis between sites. Our study supports the analysis of chlorophyll and turbidity hysteresis to infer sources and transport of suspended algal biomass.

) and eutrophication (Dodds 2006). While the aggregated behavior of phytoplankton (biomass and production) have become an integral part of river ecosystem theory and management, individual species or groups with distinct life histories may generate unique individual patterns of phytoplankton structure and function. Life history traits affect the source, transport, and fate of phytoplankton within river food webs (Reynolds 2006). Therefore, interpretations of ecosystem processes involving phytoplankton communities would likely benefit from partitioning the role of different phytoplankton groups in generating aggregated patterns of algal biomass, primary production, and respiration.
Riverine phytoplankton are typically divided into two putative groups with contrasting life histories-benthic organisms that are occasionally washed into the plankton (defined as tychoplankton), which can represent a significant portion of phytoplankton at high flows (Istvanovics et al. 2010)-and "true-plankton" species capable of actively growing while drifting. Istvanovics and Honti (2011) also hypothesized that self-sustaining plankton in rivers is predominantly meroplanktonic, that is, species adapted to spend part of their life cycle in benthic environments. Taxonomic characterization of riverine phytoplankton demonstrates the contribution of benthic and true-plankton sources to phytoplankton communities, with identifiable common taxa of Bacillariophyta, Chlorophyta, and Cyanobacteria among others (Reynolds et al. 1994;Rojo et al. 1994;Reynolds 2006). However, sorting species by life history traits is time consuming and trait overlap between groups adds uncertainty when categorizing a riverine phytoplankton community strictly by species counts.
Here, we argue that by exploring the behavior of chlorophyll concentration during stormflows (namely C-Q relationships) we can infer location and character (benthic vs. planktonic) of contributing sources to riverine phytoplankton. In the past, C-Q relationships have been widely employed to infer source areas and transportability of nutrient and sediments, but their applicability to phytoplankton remains largely unexplored (but see Istvanovics et al. 2014 andDolph et al. 2017). Our approach leverages previous interpretation of hysteresis behavior in C-Q plots for sediment and nutrients (Klein 1984;Williams 1989;Butturini et al. 2008), and relies on existing indices that quantify the rotational pattern (hysteresis index [HI]) and relative change (flushing index [FI]) in concentrations.
We expect a clockwise hysteresis over stormflow when phytoplankton sources are proximal (originate close to the point of detection and mobilize rapidly) and anticlockwise hysteresis when sources are lagged (i.e., travel time is long or mobilization delayed; Fig. 1). Hence, positive HI values indicate a phytoplankton contributing area in the nearby channel itself, including benthic (tychoplanktonic or meroplanktonic) and suspended sources that become rapidly entrained at early stages of flooding events by low flow velocities. Conversely, negative HI values are produced by events carrying lagged sources of flow-resistant benthic algae, or planktonic algae from storage-zones (e.g., backwaters) within the channel and floodplain (wetlands, ponds, side channels), that only become hydrologically connected to the channel at later stages of the hydrograph when a high threshold velocity is reached (Fig. 1). In contrast with rotational patterns of hysteresis loops, the flushing index captures the direction of the relative change in chlorophyll concentration during stormflows. We interpret that increasing phytoplankton concentrations over the rising limb (positive FI) indicate abundant benthic and/or suspended stocks of algae that can be mobilized, whereas negative FI values indicate phytoplankton dilution resulting from limited supply of algal biomass to be entrained by high flows (Fig. 1).
To test these ideas, we examined phytoplankton C-Q relationships in two locations of a river with similar discharge but contrasting hydraulic geometry due to the presence of several low-head dams ($ 1 dam per river km) in between these two locations. River channels become generally deeper and wider near dams (Stanley et al. 2002;Csiki and Rhoads 2014;Fencl et al. 2015), increasing water residence time (Zaidel et al. 2021). Our premise is that "free-flowing" river sections will favor benthic growth due to shallow waters and high light availability, leading to high phytoplankton recruitment from the benthos and thus of tychoplanktonic and meroplanktonic species. In contrast, river sections influenced by consecutive low-head dams may provide additional hydraulic storage (Reynolds and Descy 1996) favoring planktonic over benthic growth due to increased water depths and light suppression near the streambed.

Study site and field collection
We monitored high flow events at two locations in the lower part of the Brandywine Creek (Southeastern Pennsylvania), a 6 th order watershed dominated by deciduous forest (31%), pasture/ crops (31%), and low-intensity developments (24%; NLCD 2011). Our upstream site, hereafter referred to as UP Site , is located at the border of the Piedmont plateau (39.8697, À75.5934) and is adjacent to the USGS gage number 01481000 with a mean annual discharge of 9.1 AE 0.5 m À3 s À1 . Our downstream site (DW Site ; 39.76014, À75.5565) is located only 2.5-km downstream of the USGS gage number 01481500 which records a mean annual discharge of 10.2 AE 0.5 m À3 s À1 . DW Site coincides with the end of an abrupt drop in elevation characteristic of the Atlantic Seaboard Fall line. A total of 10 low-head dams still exist in the 11-km long section of the Brandywine Creek between UP Site and DW Site , compared to only two dams in the first 11 km upstream of UP Site .
At each site we installed a submersible fluorescence sensor (Cyclops 7-F; Turner Designs) calibrated for the detection of suspended Chlorophyll a concentration [Chl a]. We installed chlorophyll sensors at UP Site and DW Site on March 2019 and June 2019, respectively, and record measurements at a 5-min interval. Fluorometers were periodically cleaned and recalibrated in situ (usually every 2 weeks). We collected multiple 1 L water samples (50 at UP Site and 28 at DW Site ) to analyze acetone-extractable [Chl a] in the laboratory following standard methods (Arar and Collins 1997) in an AU-10 portable fluorometer (Turner Designs). We collected discrete samples at different stream flow levels to increase the concentration range at which in vivo (sensors) and acetone-extracted [Chl a] were compared.

Sensor data analyses
We manually delineated stormflow events by assuming that an event starts when discharge first rises and ends when discharge returns to prestorm baseflow values or when another stormflow event begins. Prior to data analysis, we processed all sensor raw data using a quality assurance and control protocol consisting of outlier removal, fouling and calibration drift correction (before-after sensor cleaning and in situ recalibration), and data gap filling by using a LOESS function to interpolate minor data gaps created during sensor maintenance and/or outlier removal. The corrected-[Chl a] values and acetone-extractable [Chla a] were highly correlated with an RMSE of 0.39 and 0.45 μg L À1 at UP site and DW site , respectively (Fig. S1).
For each event, Chl a data were paired with dissolved oxygen concentration (DO), discharge, and turbidity data from corresponding USGS sites at 15-min intervals. We estimated hysteresis and flushing indices for both Chl a and turbidity data using methods described in Lloyd et al. (2016) and Vaughan et al. (2017) developing previous methods from Butturini et al. (2008). First, concentration (C as Chl a or turbidity) and discharge (Q) values, i, for each stormflow event were normalized (x i À x min )/(x max À x min ) to allow comparisons among sites and dates. Then, we calculated the corresponding C i,norm values for Chl a and turbidity at 1% intervals of Q i,norm using a linear interpolation between the two most adjacent measurements. The hysteresis index (HI) at each discharge interval is determined by simple subtraction of corresponding interpolated values in the falling (C j,fall ) and rising (C j,rise ) limb. Stormflow HI represents the average of HI j values ( Fig. 1) and varies from À1 to 1, with negative and positive HI indicating anticlockwise and clockwise hysteresis, respectively. The magnitude of HI is proportional to the concentration difference between the rising and falling limbs (i.e., area within the loop). Then, we calculated the flushing index (FI) as indicated in Fig. 1. This index also varies between À1 and 1. Negative FI values indicate a diluting effect on concentrations during the rising limb and positive FI values indicate increasing concentrations during rising discharge, a concentration effect.We used DO and temperature data to estimate daily average rates of ecosystem metabolism (fixation and Fig. 1. The four quadrants of phytoplankton export as depicted by hysteresis (Lloyd et al. 2015) and flushing (Vaughan et al. 2017) indices. Green arrows represent predicted changes in chlorophyll concentration over time/discharge. Each quadrant represents a distinct, idealized scenario with contrasting phytoplankton sources, supply, and mobilization. Events that rapidly mobilize abundant algal biomass from in-channel sources (benthic or suspended) appear in the upper right quadrant, which contrast with events with sufficient flow magnitude to mobilize flow-resistant benthic algae and/or disconnected planktonic (e.g., backwaters, riparian wetlands, etc.) sources (upper left). Chlorophyll concentration can thus equally depict conditions dominated by benthic or true-plankton species. Events with the capacity to create phytoplankton dilution indicate a limited supply of biomass to be entrained by high flows (lower quadrants), which are most likely indicative of phytoplankton communities dominated by true-plankton species. mineralization of C) and gas exchange over the baseflow period in between high flow events using the streamMetabolizer package (Appling et al. 2018). The number of days for which we estimated river metabolic rates before each event varied from 1 to 20 d, depending on high-flow event frequency and data availability. Photosynthetically active radiation was obtained from a monitoring site 16-km away from our study sites. Water depth was estimated from discharge data using the function calc_depth() in streamMetabolizer. We used maximum-likelihood estimation to solve the model for its three parameters: gross primary production (GPP), ecosystem respiration (ER), and reaeration coefficient (K 600 ). We set the initial parameter values for GPP at 2 g O 2 m À2 d À1 , ER as À2 g O 2 m À2 d À1 , and K 600 as 2 d À1 . To assess potential equifinality problems with the model (Hall et al. 2016;Appling et al. 2018), we carefully examined ER-K 600 covariance for each model run, checked unrealistic estimations (e.g., positive ER or negative GPP), and ensured our K 600 estimations were credible based on previously published estimations (Hall and Ulseth 2020). All data used to calculate chlorophyll hysteresis indices and daily metabolism rates are available in the EDI portal https://doi.org/10. 6073/pasta/8e3fecf186d90352a05c3f63c4d793fc (Peipoch and Ensign 2022).

Results
Average daily discharge at DW Site (14 AE 4 m 3 s À1 ) was only 8% higher than at UP Site over the 2-yr study period. We obtained high-frequency Chl a data for 46 high flow events at UP Site and 24 at DW Site , with 11 of these events including simultaneous data at both locations. Event peak discharge was similar among the two sites and seasons (Table 1), varying between 9 and 183 m 3 s À1 at UP Site , and between 14 and 234 m 3 s À1 at DW Site . Average flow-weighted [Chl a] and turbidity values per high flow event were also similar between the two locations (Table 1). We found highest flow-weighted [Chl a] in spring events, followed by fall, winter, and summer regardless of location (Table 1).

Spatiotemporal variation in hysteresis behavior
Stormflow HI values of Chl a were lower on average at the UP Site , especially during the summer (Table 1; Fig. 2A). All the events at UP Site displayed positive FI values and 88% of them had a negative HI ( Fig. 2A) showing dominant concentration effects on Chl a along with counterclockwise hysteresis (lagged sources). In contrast, a slight majority (58%) of events at DW Site generated positive HI values indicating a greater role of proximal Chl a sources compared to UP Site ( Fig. 2A). Most of downstream events displayed positive FI values (concentrating effect) but 25% of them generated negative FI values (diluting effect), particularly those in which peak discharge was only 1.5 to 3 times higher than pre-stormflow levels (Fig. S2). Only 29% of the events at DW Site had negative HI-Chl combined with positive FI-Chl, a combination that occurred almost 9 out of 10 storms at UP Site . Unlike Chl a, we did not find significant differences in HI-Turbidity values between sites or among seasons ( Fig. 2B; Table 1). Counterclockwise hysteresis dominated turbidity responses at both UP Site (77% of events) and DW Site (86% of events). Values of FI-Turbidity were positive for all the recorded events at both sites (Fig. 2B), but slightly lower at DW Site in the summer (Table 1). Overall, turbidity hysteresis behavior showed prevailing concentration effects and lagged sources at both sites, which contrasts with site-specific hysteresis behavior of Chl a.

Cross-site comparison of chlorophyll and turbidity hysteresis
A comparison of Chl a and turbidity indices shows that hysteresis behavior of Chl a was akin to turbidity hysteresis at UP Site but not at DW Site ( Fig. 3A; green symbols generally lie closer to the 1 : 1 line than blue symbols). On average, HI-Turbidity changed little between sites during an individual stormflow event ( Fig. 3A; arrows connect sites measured during the same event, and most arrows point up instead of laterally). This shows that the timing of sediment entrainment was similar between sites while the proximity and contribution of chlorophyll sources differed between sites (Table 1). Similarly, FI-Chl values were generally lower than FI-Turbidity ( Fig. 3B; most values lie below the 1 : 1 line); during individual events FI-Chl usually decreased between UP Site and DW Site while FI-Turbidity stayed the same (Fig. 3B). This shows that the source of Chl a shifts from having a concentrating effect at UP Site to a more diluting effect at DW Site , even while turbidity-generating material (sediment) has a concentrating effect at both sites.

Ecosystem metabolism and chlorophyll hysteresis
Higher HI-Chl values suggest an increasing role of proximal Chl a sources that are rapidly mobilized by stormflow events (Fig. 1). To explore this causal relationship, we calculated P/R ratios from GPP and ER estimations over the baseflow period preceding stormflow events. We used net ecosystem production (NEP) and P/R ratios as indicators of the sign and magnitude of benthic algal biomass accrual due to their temporal covariance with Chl a abundance in mid-order channels (Uehlinger and Naegeli 1998). Unfortunately, our model frequently converged to high K 600 rates (>20 d À1 ) and/or positive ER values when using data from DW Site, and this prevented us from estimating daily metabolic rates at DW Site -possibly due  to the effects of dams on reaeration. At UP Site , daily GPP averaged 2.3 g O 2 m À2 d À1 and varied from 0.33 to 5.28 g O 2 m À2 d À1 , mean ER was À3.03 g O 2 m À2 d À1 and K 600 averaged 3.45 d À1 throughout the study period. NEP rates and P/R ratios prior to high flow events had significant effects on chlorophyll concentrations and hysteresis behavior. We found that average daily P/R ratios prior to stormflow events were positively correlated with the event's flow-weighted chlorophyll concentrations (Fig. 4A). Similarly, both P/R and NEP also predicted event HI-Chl values at the UP Site (Fig. 4B). This indicates that pre-stormflow periods with autotrophic conditions (i.e., more positive NEP and P/R) increased the contribution of proximal Chl a sources into Chl a hysteresis behavior. At Up Site , this occurred predominantly during spring and fall periods of sparse canopy and maximum irradiance (Fig. 4). We found no relationship between P/R ratios and HI-Turbidity (R 2 = 0.02; b = 0.09; p-value = 0.72).

Discussion
Hysteresis behavior of phytoplankton varied along a river segment reflecting reach-scale hydrogeomorphic characteristics (i.e., abundance of low-head dams) and local sources of primary production. Concentrating effects on suspended chlorophyll at the free-flowing section suggest that an abundant supply of tychoplankton and meroplankton species were being entrained from the benthos after a velocity threshold was exceeded and that induced a delay in hysteresis. When NEP preceded high flow events in the spring and fall, chlorophyll lag decreased and export increased, likely due to a higher benthic biomass contribution of growth forms (e.g., filamentous algae) that are less flowresistant but characteristic of autotrophic riverine conditions (Power 1990). Conversely, chlorophyll was less concentrated or even diluted by storm flows downstream of a series of low head dams. This was most notable during low magnitude floods capable of minimum algal re-suspension but significant dilution of true-plankton concentrations. We presume that the additional hydraulic storage, water depth, and resultant light suppression caused by the low-head dams reduced tychoplankton biomass and favored true-plankton contribution. The contrasting effects of storm flows on chlorophyll concentration at each location compared to their highly similar concentrating effects on turbidity (i.e., sediment) also suggests the existence of Chl-depleted but sediment-rich benthic sources near the dams. However, a massbalance estimation of the necessary suspended chlorophyll concentration in storage-zones (Reynolds and Descy 1996) suggests that true-plankton growth in the lower section of the watershed was not sufficient to generate the observed Chl a export. Even if we assume that storage-zones (i.e., low head dams) occupied 40% of the water volume between our two sampling locations, the Chl a concentration in them would need to be 5-300 times greater than baseflow Chl a concentrations at DW Site to account for the observed Chl a export during stormflows. This is unlikely. Instead, we believe the phytoplankton community in the Brandywine Creek is mostly sustained by benthic recruitment (Istvanovics and Honti 2011) with a moderate influence of trueplankton growth in the lower section of the watershed.
Our case study focused on two contrasting locations with obvious geomorphic differences supports the interpretation of chlorophyll hysteresis patterns to discriminate between benthic and suspended algae. We acknowledge that further development of this technique should focus on how thresholds in shear stress, hydrologic connectivity of off-channel areas, and biological heterogeneity and transport at the reach-scale influence interpretation of hysteresis. Nevertheless, phytoplankton groups with contrasting life histories (Istvanovics and Honti 2011), taxonomy (Reynolds et al. 1994;Rojo et al. 1994), and metabolic activity (Dokulil 2014) must contribute Seasonal variation and NEP are also illustrated in each plot. Dashed line indicates the frontier between heterotrophic (P/R < 1) and autotrophic (P/ R > 1) river metabolism.
differently to aggregated patterns of community structure and ecosystem function. For example, riverine production dominated by meroplankton might reset more strongly following floods than one dominated by a true-plankton source upstream, information critical for interpreting time series trends in river metabolism. In fact, our results showed that Chl a hysteresis not only reflects the distinct contributing algal sources but also the availability and transport distances, or time lags, of autochthonous C within the ecosystem boundaries. These differences could affect how phytoplankton (organic carbon) fluxes are attributed to source areas and the factors driving algae growth rates, and hence affect regulatory measures aimed at controlling algal growth.
Despite significant advances in O 2 -based models as a proxy for river metabolism (Appling et al. 2018;Hall and Ulseth 2020), the underlying factors driving autochthonous production and fate across watersheds of different size are not yet fully understood (Thorp et al. 2006;Bernhardt et al. 2018). Coupling ecosystem-level metabolism measurements to Chl a hysteresis across rivers of contrasting size, hydrogeomorphic complexity, and/or human influence may help fill existing knowledge gaps on river autotrophy and upstream-downstream energetic linkages. Indeed feedbacks between the river's physicochemical state and metabolic processes during high flow events have recently been documented (O'Donnell and Hotchkiss 2019). We contend that better understanding of aggregated ecological function, namely nonstationary attributes (GPP or ER rates), might also depend on greater characterization of ecosystem state variables such as phytoplankton composition and benthic vs. plankton biomass pools. For instance, our results suggest that, depending on water-column depth and light availability, the metabolic activity of benthic autotrophs may be significantly reduced if benthic biofilms are being subsidized by true-plankton deposition. In that case, correlational approaches (McTammany et al. 2007) or temporally dynamic models (Segatto et al. 2020) linking metabolic rates and benthic autotrophic biomass are not likely to capture the true ecological state of the benthic algal community in mid-to-large river sections.
Indeed, the riverine productivity model (Thorp and Delong 1994;Thorp and Delong 2002) argues for the disproportionate importance of autochthonous C in heterotrophic rivers due to localized but intense autotrophic activity. This is particularly applicable to mid-and high-order rivers where phytoplankton communities can develop in the main channel or in connected slackwaters (Thorp et al. 1998;Lewis et al. 2001), but evidence shows that autochthonous C can be a primary energy source for secondary production in an array of stream orders, including oligotrophic, low-order systems draining well-forested catchments (Hill et al. 2001;Mulholland et al. 2006;Roberts and Mulholland 2007). Our results pose a promising approach with minimal methodological limitations (appropriate for many stream sizes) to further investigate spatial and temporal variation in source and transport of autotrophic biomass, although further research linking chlorophyll hysteresis patterns to algal taxonomy, cell size, or Chl/C ratios is surely necessary.