Transport of biodeposits and benthic footprint around an oyster farm, Damariscotta Estuary, Maine

Study area

The Damariscotta is a narrow, partially mixed mesotidal estuary in the midcoast of Maine. The study was conducted at Perkins Point, approximately 20 km from the mouth of the Damariscotta River (44 00.016′N, 69 32.650′W), where a large oyster aquaculture operation is located on a 24-acre (∼97,000 m²) lease (Mook Sea Farm, Fig. 1). Verbal permission to use the site for the study was obtained from Bill Mook, owner/operator of Mook Sea Farm. Eastern oysters (Crassostrea virginica) are grown at the farm in mesh bags held in rigid, rectangular floating OysterGro™ cages that are strung together in rows. The farm sits above a sloping shoal with a mean water depth ranging from ∼2–4 m. We assessed our estimate of the footprint by comparing a site directly underneath the oyster cages (“Farm” site) in the southern portion of the farm with one downstream, ∼90 m south of the edge of the farm lease area (“Away” site). The Away site location was constrained by the geomorphology of the river; because the river bends and narrows downstream (Fig. 1), selecting a site farther than ∼90 m would risk it being in an area too hydrodynamically different from the Farm site for comparison. Instruments were deployed on moorings at both sites to measure hydrographic properties, with sediment traps to quantify sedimentation. Sediment cores were also collected at each site for evaluation of sediment erodibility, geochemical characteristics, biogeochemical fluxes, and macrofaunal community structure. While diving to collect the cores for this study, we observed benthic microalgal mats in distinct patches of varying sizes (Fig. 2), though we did not discriminate between areas with and without mats in our sampling.

Flow and deposition

To determine the horizontal distance particles would be transported before deposition (Q1), flow velocity near the bottom was measured over a tidal cycle with an Acoustic Doppler Velocimeter (ADV) and combined with laboratory measurements of particle settling velocities. Here we define the footprint based on the settling velocity of the mean size of biodeposits and the maximum distance transported over the tidal cycle during neap tide. Because mean water column height varied from ∼2–4 m across the farm, calculations were done with the maximum and minimum depths to provide a range. Thus, biodeposit particles smaller than the mean would travel beyond our defined footprint under peak flow velocities, but the footprint would encompass most of the tidal range in deposition. To demonstrate the effect of particle size on the footprint, we also calculated transport distance for the 5th and 95th percentile of particle sizes.

The ADV was horizontally mounted to a metal tripod mooring deployed at the Farm site, measuring velocity (8 Hz) 37 cm above the bottom. An Onset HOBO dissolved oxygen (DO) logger was also mounted to the mooring (50 cm above the bottom) to measure near-bottom DO to determine whether low oxygen occurred that might affect biogeochemical fluxes or infauna. A cinder block mooring was deployed at the Away site. Temperature and salinity measurements were collected near-surface and at bottom at both sites to measure stratification. Water level at the farm was measured by the pressure sensor on the ADV. Moorings were deployed for ∼23 h on August 2–3, capturing nearly two tidal cycles during a neap tidal period. Because the experiment was conducted over a short time frame, hydrographic data for the entire summer season was obtained from a Land/Ocean Biogeochemical Observatory (LOBO) surface mooring ∼250 m southeast of the farm (part of the Sustainable Ecological Aquaculture Buoy Network; http://maine.loboviz.com).

Sediment traps were attached to the bottom moorings at both sites to measure sedimentation rates. The traps were made of 5.08 cm inner diameter PVC cylinders, 15 cm in length, and capped on one end, following Gardner (1980). Traps on the bottom mooring at the Farm site were diver-deployed and intended to measure a combination of biodeposits directly from the farm as well as ambient deposition. The Away site traps were attached to the cinder block mooring that was lowered to the bottom from the boat and were intended to measure the ambient particle deposition without influence from the farm. The mooring was repositioned once it was on the bottom which may have created a sediment plume that falsely inflated sedimentation data at that site, so those data should be interpreted with caution. Additional traps were attached directly to the underside of the oyster cages at the Farm site for a direct estimate of in-situ particle production by the oysters without under-farm flow transporting the particles. Comparison of these data to deposition at the bottom mooring took into account the area of the farm occupied by gear (∼40% - obtained from the farm’s lease application). Six traps were deployed during two 24-hr deployment periods separated by 5 days, during which the tidal range varied by ∼0.5 m. All traps were deployed and retrieved at high tide. At least one trap was lost from each location during each deployment. During the first period, traps were deployed and recovered from the Farm site (n = 4) and the cages (n = 5). During the second period, sediment traps were deployed and recovered from the Away site (n = 5) and the cages (n = 5). After collection, sediments were dried for 24 h at 60 °C and weighed.

To measure biodeposit settling velocities, biodeposits were collected from adult, harvest-size oysters held in a flow-through seawater system. Biodeposits were spread in a petri dish to prevent particle aggregation, maintaining the size and cylindrical shape of freshly released biodeposits observed in both the lab and field. Deposits were then pipetted into the middle of a 37.8 L aquarium of seawater, at least 12 cm away from any wall to prevent wall effects. Videos of their descent were recorded using a Nikon D5300 DSLR camera equipped with an AF-S Micro NIKKOR 105 macro lens. It was not possible to distinguish between feces and pseudofeces from the observations, and we did not describe the composition of the biodeposit particles. The descents of 100 individual particles were recorded. Each particle’s dimensions (length and diameter of cylindrical particle) were measured in ImageJ and the angle (normal to vertical) and rate of descent tracked with the ImageJ MtrackJ plugin (Schneider, Rasband & Eliceiri, 2012; Meijering, Dzyubachyk & Smal, 2012). Settling velocity depends on density and particle surface area. The length and width of each particle measured from the videos were used to calculate surface area using the geometric formula for cylinder surface area. We were unable to measure the densities of individual pellets, but density was calculated using a form of Stokes’ settling equation specifically modified for small cylinders (Komar, 1980), (1)${\displaystyle V_t=0.079\frac{1}{\mu }\left({\rho }_p-{\rho }_w\right)g{L}^2{\left(\frac{L}{D}\right)}^{-1.664}}$ where V_t is settling velocity, µ  is the viscosity of the water, ρ_p and ρ_w are the densities of the particle and water respectively, g is acceleration due to gravity, L is particle length, and D is particle diameter. Particle surface area and angle of descent were compared to settling velocity using a multiple regression. Density was not included in this model, as it was calculated from the settling velocity.

An estimate of biodeposit tidal transport (L_f) was calculated as: (2)${\displaystyle L_f=\frac{u_s\left(t\right)h\left(t\right)}{w_s}}$ using the bottom horizontal velocity (u_s) measured with the ADV and the mean particle settling velocity (w_s) calculated from the lab observations. The height of the water column (h(t)), or vertical distance traveled, and bottom horizontal velocity (u_s(t)) changed depending on the time (t) in the tidal cycle. To provide an upper and lower bound on the potential biodeposit tidal transport, calculations used the deepest and shallowest mean water depths at the farm (2 and 4 m). It should be noted that this calculation uses velocity measured with the ADV which is assumed to be in the bottom boundary layer, so the transport distances should be viewed as conservative estimates.

Erosion

To determine whether settled biodeposits were likely to be resuspended (Q2), we calculated bottom shear velocity throughout the flood/ebb tidal cycle and conducted laboratory erosion experiments across a similar range of bottom shear velocities.

We used the flow velocities measured by the ADV to calculate bottom shear velocity throughout the tidal cycle via three indirect methods (logarithmic profile, covariance and turbulent kinetic energy methods) as described by Kim et al. (2000) for high sampling rates (>5 Hz). Multiple methods are commonly used and compared, as methods can be biased by local conditions such as waves, stratification, and bedforms (e.g., Sherwood, Lacy & Voulgaris, 2006; Pieterse et al., 2015). The logarithmic profile method assumes a logarithmic velocity profile using the von Karman-Prandtl equation. The reliability of the estimates from this method is severely limited by having only a single point velocity measurement. However, single point techniques have been developed (Sherwood, Lacy & Voulgaris, 2006) and Kim et al.’s (2000) methodology was followed. The covariance method requires the velocity measurements to be in the constant stress layer of the bottom boundary (a reasonable assumption for the instrument depth and low water column stratification, buoyancy frequency ∼0.004) and uses turbulent fluctuations in the along-estuary direction and vertical flow velocities to determine the bottom stress. The turbulent kinetic energy method assumes a linear relationship between turbulent kinetic energy and bottom stress. Bottom shear stress (τ_b) from both of these methods is then linked to shear velocity as ${\tau }_b=\rho u_{\ast }^2$ where ρ is water density and u_∗ is shear velocity. As seen in other studies, all three methods produced results with similar patterns, but peak magnitudes showed considerable differences (e.g., (Kim et al., 2000; Sherwood, Lacy & Voulgaris, 2006; Pieterse et al., 2015)). Because the covariance method is less affected by local conditions (e.g., stratification, bedforms; (Sherwood, Lacy & Voulgaris, 2006)) and the relationship between the along-estuary velocity and shear velocity was well represented by a quadratic fit ($u_{\ast }=0.27096u_s^2+0.00089u_s+0.00466,R=0.93$), we used that estimate for our further analyses. The effects of waves on bed stress were not considered because of the short fetch of the estuary.

We conducted laboratory erosion experiments using a custom-built Gust erosion chamber to generate near-uniform bed shear (Gust & Muller, 1997; Thomsen & Gust, 2000; Tengberg et al., 2004; U-GEMS Manual, Green Eyes LLC, 2015). Triplicate 10 cm diameter sediment cores were collected by divers from the Farm and Away sites and placed in a large holding tank with flowing seawater from the Damariscotta River. All erosion tests were performed within 60 h of collection. Each core was capped with the Gust chamber, and a rotating disc within the chamber generated increasing levels of shear velocity (0.30 cm s⁻¹, 0.95 cm s⁻¹, 1.34 cm s⁻¹, 2.01 cm s⁻¹, 2.32 cm s⁻¹). Each shear velocity level was maintained for 20 min before increasing disc rotation and flow rates. For each 20 min shear velocity level, water and eroded material were removed by a pump at the center of the disc at 10 cm above the sediment surface and filtered through 47 mm Whatman GF/F filters (1.5 µm pore size). Effluent was replaced with the same seawater supplied to the core holding tank, so suspended mass in the pre-erosion core overlying water and the replacement water was likely similar. The lowest shear velocity, 0.30 cm s⁻¹, was used as a flushing step and was not filtered for suspended sediment analysis (no resuspension was observed at this low shear). Filters were then dried at 60 °C for 24 h and weighed to obtain eroded sediment mass for each shear velocity level. One Away core was not used in the analysis due to accidental introduction of a large nereid worm from an unrelated experiment in the holding tank. Suspended sediment concentration, C_s (kg m⁻³), for each site at each velocity level was calculated from the dry mass (kg) of filtered sediment divided by the volume (m⁻³) of water filtered. C_s was converted to eroded mass per area (E; kg m⁻²): (3)${\displaystyle E=\frac{C_s{V}_c}{A_c}}$ where V_c is chamber volume (7.9 ×10⁻⁴ m³), and A_c is sediment surface area within the core (7.9 × 10⁻³m²). Due to small and unequal sample sizes, statistical comparisons between sites and stress levels were not performed.

To generate specific shear velocities, cap rotation and pumping rate were set using calibration equations from the University of Maryland Center of Environmental Science Gust Erosion Microcosm System (U-GEMS) Manual (Green Eyes LLC, 2015): (4)${\displaystyle u_{15}^{\ast }=0.0318n^{0.763}}$ (5)${\displaystyle Q=-28.31u_{15}^{\ast 2}+170.2u_{15}^{\ast }-23.85}$ where u^∗₁₅ is shear velocity at 15 °C (cm s⁻¹), n is rotations per minute, and Q is pumping rate (mL min⁻¹). Shear velocity at 15 °C was converted to shear velocity at the water temperature measured during the erosion tests (22 °C) as: (6)${\displaystyle u_{15}^{\ast }=u_{22}^{\ast }\left[1+0.006\left(22-15\right)\right]}$ where u^∗₂₂ is shear velocity at 22 °C (cm s⁻¹) (U-GEMS Manual, Green Eyes LLC, 2015). For those more familiar with stress as a measure of shear, the corresponding applied bottom shear stresses (τ_b; Pa) generated by the Gust chamber were 0.009 Pa, 0.09 Pa, 0.18 Pa, 0.41 Pa, and 0.55 Pa. These were calculated from shear velocities (u^∗₂₂; m s⁻¹) as: (7)${\displaystyle {\tau }_b=\rho u_{22}^{{\ast }^2}}$ (U-GEMS Manual, Green Eyes LLC, 2015). A ρ of 1,021 kg m⁻³ was used based on an average temperature of 22 °C and salinity of 31 psu.

This field sampling occurred during neap tide, and also coincided with the highest surface water temperatures of the year when thermal stratification likely peaked. To gain a broader sense of the potential for sediment resuspension during other times in the tidal cycle, near-surface turbidity and velocity data from the LOBO surface mooring east of the farm were analyzed for longer-term patterns. While turbidity is not a direct measure of sediment resuspension, water clarity is typically related to sediment concentrations. Potential relationships among precipitation, water level, current velocity and turbidity were investigated using both the instantaneous data and 40-hour low-pass filtered signals. Near-surface velocity data were also used to estimate bottom shear velocities in higher-flow times in the spring-neap cycle. Near-surface velocity was converted to a proxy for near-bottom velocity based on a scaling factor determined from the relationship between the LOBO buoy and the ADV-measured velocities during the overlap period in August. This near-bottom velocity proxy was used to estimate shear velocity using a quadratic relationship between the longitudinal velocity and direct calculations of the shear velocity (using the covariance method). This fitted relationship was extrapolated to predict shear velocity expected from the highest flow in the spring-neap tidal cycle.

Geochemical characterization

Sediments at both sites were evaluated for grain size distribution, water content, percent organic content, and concentrations of C and N. One diver-collected 10-cm diam core from each site was vertically sectioned. Water content was calculated as the difference between wet and dry masses after sediment was dried at 60 °C for 24 h. Percent organic content was calculated as the difference between dry and ash masses after burning sediment at 500 °C for 4 h. Burned sediment was then soaked in a 1% sodium metaphosphate dispersant solution for 24 h to rehydrate and disaggregate sediment prior to grain size analysis. Grain size distribution was determined using a Malvern Mastersizer 3000 particle analyzer (Malvern Panalytical, Malvern, UK). Data were analyzed using Gradistat (Kenneth Pye Associates, Ltd., Berkshire, UK). Preparation of C and N samples included an acidification step to remove inorganic carbon. Concentrated HCl was added to dried sediment samples and heated at 80 °C until the fluid evaporated as described by EPA (2009). Acidified and non-acidified samples were analyzed for C:N to determine the proportion of inorganic carbon. C and N were analyzed on a CHNSO elemental combustion system analyzer (mod. ECS 4010; Costech Analytical Technologies, Valencia, CA). C:N is reported as a molar ratio (mmol g⁻¹ sediment).

Biogeochemical fluxes

To determine whether the “Away” site exhibited similar impacts of oyster biodeposits on biogeochemical fluxes as the site directly under the farm, sediment core incubations were performed in the dark to measure nutrient, DIC and oxygen fluxes and sediment surface samples were collected for chlorophyll a (chl-a) analysis. In an additional treatment, “Amended” cores, we added oyster feces and pseudofeces onto the top of two cores from the “Away” site. Because of the lack of available data for this area on typical respiration rates, the Amended treatment was included to determine whether Away site sediments were being enriched with biodeposits; if sediments from the Away site were receiving large influx of organic matter from the farm, i.e., were saturated, additional organic matter would be expected to have minimal effect on nutrient and oxygen fluxes. The usefulness of the Amended cores therefore depended on “background” (i.e., unenriched) respiration rates being relatively low.

Incubations of sediment cores were performed to measure sediment oxygen flux and nutrient (NH₄⁺ and PO₄³⁻) and dissolved inorganic carbon (DIC) fluxes. Only dark incubations were performed so that maximum respiration rates could be measured without the influence of microalgal photosynthesis. Our measurements are therefore not representative of daytime in situ fluxes. A total of 12 15.2-cm diam sediment cores were taken by divers during two collections (2 and 7 August 2018). During each collection, two cores were taken from the Farm site and four were taken from the Away site. Care was taken during collection that each core had at least 10 cm of sediment depth and at least four cm of overlying water. The cores were brought back to the lab and kept in tanks with flowing water from the Damariscotta River until use in the incubations (<36 h post-collection). Immediately before the start of the incubations, eight mL wet volume of oyster feces and pseudofeces were pipetted onto the top of two cores collected from the Away site and designated the Amended cores. Between the two samplings, this resulted in a total of four replicates each of Farm, Away and Amended cores.

In incubations, the cores were each submerged in separate 5-gallon buckets of seawater and watertight caps affixed to the tops (after methods and Fig 2 in Dorgan et al., 2020). The caps allowed for six cm of overlying water and were each equipped with a stir-bar to gently circulate the overlying water, inflow and outflow taps for water sampling, and an Atlas Scientific DO probe (Atlas Scientific LLC, Long Island City, NY) that measured oxygen concentrations in the overlying water of each core every ∼12 s. The buckets containing the cores were covered with foil to prevent oxygen production by photosynthesis from corrupting oxygen consumption rate measurements. The incubations of cores collected on 2 August provided simultaneous oxygen, nutrient and DIC fluxes. Nutrient and DIC samples were gathered during the incubation of the 7 August cores, however the oxygen measurement instrument failed to record during the initial incubation, so a second incubation was done on the same cores the following day to collect data to calculate oxygen flux. Sediment samples for chl-a were taken from the top ∼1cm of the cores following the second incubation and frozen for later analysis. The samples were prepared using methods adapted from Welschmeyer (1994) and analyzed on a Turner Trilogy 7200 using an HCl acidification step to measure chl-a and phaeopigment concentrations.

Nutrient and DIC samples were siphoned from overlying core water at three timepoints using methods described by Lehrter et al. (2012). Samples were taken at the start (0 h elapsed), middle (2 h elapsed) and end (4 h elapsed) of the incubations. DIC samples were kept dark and cold in 20 mL glass vials until processed (<24 h post-collection) on a TOC Carbon analyzer. Nutrient samples were filtered through 25 mm Whatman GF/F filters (0.7 um pore size) and frozen in the dark until processed for NH₄⁺ and PO₄³⁻ at the Dauphin Island Sea Lab. NH₄⁺ was measured fluorometrically with modifications as described by United States Enviromental Protection Agency (2012) and Holmes et al. (1999). PO₄³⁻ was measured on a spectrophotometer with slight modifications from methods described in Grasshoff, Ehrhardt & Kremling (1983). Rates of change in overlying water concentration for oxygen, DIC and each nutrient were calculated from the linear regression of the respective analyte concentrations versus time. These rates of change (mmol min⁻¹) were multiplied by the height of overlying water to obtain flux rates (mmol m⁻² d⁻¹). We were unable to obtain oxygen flux values for the Amended cores because they lacked a consistent oxygen trend. Some cores produced several successive oxygen slopes within the same range of oxygen concentrations because low-oxygen water was replaced with high oxygen water during the sampling for nutrients. To calculate a single oxygen flux for each of these cores, time-weighted average slopes were calculated. Shapiro–Wilk and Kruskal–Wallis tests were performed on nutrient and oxygen flux data to determine if there were differences in fluxes between sites (R Core Team, 2016).

Macrofaunal community structure

After data collection, all incubation cores were sieved (500 µm) and the remaining material preserved in 95% ethanol and rose Bengal stain. Macrofauna were picked out of each sample, and individuals identified to family level. Statistical comparison of abundances between the two sites was not performed due to small and uneven sample size. Data from replicate cores were analyzed for Shannon diversity (H’) and Pielou’s evenness (J’), averaged by site, then Student’s t-tests were performed on the diversity and evenness metrics to compare the two sites. An NMDS analysis was performed on square-root transformed data (to decrease weighting of highly abundant groups) to compare community structure between sites (PRIMER v7, PRIMER-E Ltd, Plymouth).

Flow and deposition

The velocity data were strongly tidal during the deployment, with a clear semi-diurnal cycle. The water level had a 2.86 m range, with a mean water depth of 3.4 m at the away site and 3.5 m at the farm site. The ADV measured near-bottom velocities ranging ± 20 cm s⁻¹ (Fig. 3A). While the magnitude of the flood and ebb currents were similar, there was temporal asymmetry in duration, with flood tide being ∼1.9 times longer than ebb. Bottom dissolved oxygen remained above 8 mg L⁻¹ throughout the sampling period (data not plotted).

The rate of deposition calculated from sediment traps attached to the bottom mooring beneath the farm was (490 ± 36 g m⁻² d⁻¹, n = 4). The cage-area normalized rates of deposition calculated from the traps attached to the cages were (163 ± 86 g m⁻² d⁻¹, n = 5) for the first deployment and (230 ± 132 g m⁻² d⁻¹, n = 5) for the second. These rates were calculated by multiplying the deposition rate calculated from the beneath-cage traps by the percent area of the farm occupied by gear (40%), following the methods of Testa et al. (2015). The highest deposition rate was at the Away site (892 ± 87 g m⁻² d⁻¹, n = 5), where deployment and repositioning of the cinder block mooring may have resuspended sediments.

The biodeposits were cylindrical in shape, with an average length of 0.64 ± 0.14 mm and average diameter of 0.11 ± 0.02 mm. Biodeposit surface area ranged from 0.12 mm² to 0.60 mm² with an average surface area of 0.25 ± 0.09 mm². For each one mm² increase in surface area, we observed a 0.57 (± 0.53; ± 95% CI) cm s⁻¹ increase in settling velocity (p = 0.03; r² = 0.03) (Fig. 4). The angle of descent did not significantly affect the settling velocity of biodeposits (p = 0.72). Using the modified Stokes equation, mean pellet density was calculated to be 1.42 ± 0.18 g cm⁻³with a range of 1.13 to 2.26 g cm⁻³.

To simplify calculations of horizontal particle transport, settling velocity was averaged across all particle sizes as 0.73 ± 0.23 cm s⁻¹. At this velocity, a particle would take ∼5-10 min from release to deposition across the range of typical site depths (2–4 m). Oyster biodeposits had a mean tidal transport range of 72 m upstream at peak flood and 62 m downstream at peak ebb in the shallow (2m depth) section of the farm (Fig. 3B) and 118 m both upstream and downstream in the deep (4m depth) section of the farm (Fig. 3C).

Erosion

Bottom shear velocities generated by the Gust chamber were in the range of those calculated from the bottom velocities measured by the ADV (Fig. 5A). Shear velocity from the covariance estimate ranged from 0.27 cm s⁻¹ during slack tide to 1.7 cm s⁻¹ at max ebb and flood. The lowest shear velocity tested in the erodibility experiments (0.95 cm s⁻¹) was within this range, and the higher shear velocities tested extended above the maximum field observations, with values up to 2.32 cm s⁻¹ (Fig. 5A, Fig. 6). There appeared to be more eroded mass at a given shear velocity under the farm than away and average total cumulative eroded mass under the farm (6.82 ± 2.52 g m⁻²) was over twice as high as that away from the farm (3.22 ± 0.34 g m⁻²). Additionally, eroded mass appeared to increase more with increasing shear velocity at the Farm than at the Away site, but only at the highest shear velocity did error bars not overlap between the two sites (Fig. 6).

There was visible evidence of erosion that was not quantified by the erosion chamber. Cohesive microbial mats were observed on the sediment surface in most erosion cores, and much of the eroded material was composed of larger particles or aggregates that were not entrained high enough in the erosion chamber to be taken up in the effluent and quantified by filtration. During the 2.01 cm s⁻¹ period of one Farm replicate, dislodged aggregates formed a mound in the center of the core (Fig. 5B). Partial mat failure was observed in the 2.32 cm s⁻¹ period of one Away replicate (Fig. 5C).

Longer-term data from the nearby LOBO buoy show that surface velocities under spring tidal conditions were nearly twice as high as during the period of the field sampling (Fig. 7). The scaling factor used to estimate bottom velocity from the relationship between the surface LOBO buoy and bottom ADV was 1.40 and showed a good fit (R = 0.95). Using the scaling factor, bottom velocities during spring tides were estimated as ∼30 cm s⁻¹ during ebb and ∼25 cm s⁻¹ during flood, with the maximum at ∼40 cm s⁻¹. Neap tide bottom velocities were ∼20 cm s⁻¹ during both ebb and flood, corresponding well with the bottom velocities measured by the ADV. During the spring tides, the calculated shear velocity approached three cm s⁻¹, higher than the highest shear velocities tested in the erosion tests or encountered during the experiment (Fig. 5A). The long-term data also show that peak velocities were often associated with high turbidity events (Fig. 7). It should be noted that some of the high turbidity events were also associated with large precipitation events (triangles in Fig. 7), thus were likely driven by sediment runoff rather than resuspension.

Geochemical characterization and biogeochemical fluxes

Sediments at both sites were poorly sorted coarse silt to fine sand (Fig. 8A, Table S1) (Folk & Ward, 1957). Organic content and water content were higher in muddier sediments as expected, and muddiness increased with depth at the Away site but varied with depth at the Farm site (Fig. 8, Table S1). The C:N (mmol g⁻¹ sediment) was 11.2 at the Away site and 11.48 at the Farm site (Table S1). Surface sediment chl-a measurements from Farm sediments were 14.3 and 11.8 µg chl-a g⁻¹ dry sediment and Away sediments were 7.3 and 22.5 µg chl-a g⁻¹ dry sediment. Phaeopigment concentrations were lower than chl-a in both samples from both the Farm site (5.88 and 6.80 µg g⁻¹ dry sediment) and Away site (11.05 and 8.74 µg g⁻¹ dry sediment).

Replication was not high enough to test if measurements of sediment fluxes differed between the two incubations, so data from both incubations were pooled by treatment (Away, Farm, and Away Amended, n = 4 for each treatment, Table S2). Sediment oxygen flux did not differ between the Farm and the Away sites, with negative values indicating flux into the sediment (Mann–Whitney Wilcoxon p = 0.69) (Fig. 9A). DIC fluxes out of the sediment did not differ among treatments (Kruskall–Wallis p = 0.66) but were highly variable, especially among the Away cores (Fig. 9B). NH₄⁺ and PO₄³⁻ fluxes were similar among all treatments (Kruskall–Wallis p = 0.79 and 0.31, respectively) (Figs. 9C, 9D). There was no relationship between fluxes of O₂ and DIC (p = 0.36), NH₄⁺ (p = 0.44) or PO₄³⁻ ( p = 0.51) (Figs. 9C, 9E, 9G).

A large nereid worm was found in a core taken from directly underneath the farm. This core also had the highest oxygen flux of all cores measured (−72.5 mmol m⁻² d⁻¹). A large nephtyid worm was found in an Away core that had the next highest oxygen flux measured (−60.6 mmol m⁻² d⁻¹). Removing these two cores from the sediment flux analyses did not result in significant differences among treatments (all p > 0.1) or significant relationships between nutrients and oxygen flux (all p > 0.45).

Macrofaunal community structure

Twelve families of macrofauna in 3 phyla were found among all samples collected (Fig. 10A, Table S3). Overall, macrofauna present at both sites were small-bodied and mainly consisted of burrowing, suspension-feeding bivalves (e.g., Myidae, Mactridae), surface deposit-feeding polychaetes (e.g., Spionidae, Flabelligeridae) and the burrowing, carnivorous polychaete, Nephtyidae. Abundance and taxa richness were highly variable at both sites. Family richness did not appear to differ between the sites, with 3.25 ± 1.89 families (mean ± st. dev.) at the Farm site and 4.0 ± 2.0 families at the Away site. Abundance appeared to be higher at the Away site, with 3,411 ± 3,902 individuals m⁻², compared to 459 ± 397 individuals m⁻² at the Farm site, but variability was high with the abundance of most families having a greater standard deviation than the mean (Fig. 10A, Table S3). Abundance was low in all macrofaunal groups except for small nephtyid polychaetes which dominated at the Away site, making up 88% of the abundance. Shannon diversity also appeared similar but variable, with H’ = 0.85 ± 0.22 at the Farm site and H’ = 0.71 ± 0.60 at the Away site (Student’s t-test, p = 0.59). Pielou’s evenness (J’) at the Farm site was 0.86 ±  0.18 and at the Away site was 0.54 ± 0.36 and was not significantly different between sites (Student’s t-test, p = 0.08). There was no distinct partitioning in community structure between Farm and Away samples, although several of the Away site samples that shared relatively high abundance (> 6,000 m⁻²) of nephtyids were > 60% similar (bold circle in Fig. 10B).

Farm footprint

The near-bottom velocities at the study site during the study period and lab measurements of biodeposit settlement rates were used to calculate a rough “footprint” of the impact of the farm on the sediment. An average-sized biodeposit released in the deeper section of the farm likely settled at a maximum of 118 m up or downstream of the farm at peak tidal flows, and far less if it was released in the shallower part of the farm (Fig. 3). At peak flow the smallest biodeposits in the deep section of the farm could have traveled up to 200 m from the farm before deposition, however transport of particles in the shallow section was more constrained (Fig. 3B), indicating that most particles would be settling within a far smaller range. Based on our footprint estimation, our selected Away site at ∼90 m downstream falls within the footprint of the farm and would therefore be expected to regularly but not continually receive oyster biodeposits.

There are few recorded measurements of oyster biodeposition settling velocities, however our mean rate of 0.73 cm s⁻¹ used in the “footprint” calculation is similar to previous measurements of ∼0.8cm s⁻¹ (Haven & Morales-Alamo, 1972) and, unlike previous investigations, was measured using biodeposits from the system of study. Our estimate of settling velocity is conservative, as it assumes biodeposits are sinking unaggregated. In situ, biodeposits may aggregate together, increasing settling velocity and reducing the farm’s “footprint”.

The weak correlation between settling velocity and surface area (Fig. 4) in addition to the non-zero intercept (we would expect particles with a surface area of 0 mm² to have a settling velocity of 0 cm s⁻¹) illustrates that either the relationship between surface area and settling velocity is nonlinear near zero, or that other parameters such as biodeposit density vary in a non-random way. We hypothesize that the scatter of data and the related poor fit likely results from differences in density among biodeposits, specifically that smaller biodeposits settle relatively faster because they are more dense than larger biodeposits. Oysters separate labile organic matter from refractory material and nonorganic material like sand and excrete these materials in feces and psuedofeces, respectively (Haven & Morales-Alamo, 1966). Biodeposits at the lower end of our range of densities (1.14 g cm⁻³) are likely composed primarily of organic matter, while those at the upper end (2.26 g cm⁻³) likely consist primarily of nonorganic minerals. Presumably, the relative makeup of these deposits depends on the composition of the suspended material in the water column at the time of filtration, affecting the densities of both feces and pseudofeces. We were unable to differentiate between feces and pseudofeces in this study, and did not measure the density of individual pellets, preventing stronger correlations between biodeposit characteristics and settling velocity.

The biodeposits of different species of bivalves also appear to have different settling velocities, and larger bivalves tend to produce more biodeposits than smaller ones (Haven & Morales-Alamo, 1972). It is therefore important to consider the species and the size of bivalves when estimating the area of impact. For example, Callier et al. (2006) found that settling velocities for marine mussel biodeposits ranged from ∼0.3–1.8 cm s⁻¹ and varied based on the size of the individual. These higher settling velocities would result in an area of impact smaller than that calculated here based on oyster biodeposits, which may lead to greater rates of accumulation, particularly if individual bivalves were larger and produced more biodeposits. We did not account for individual size in this study, which could have been a source of variability in both the mean settling velocity and mean deposition rate.

Cage-area normalized rates of deposition ranged from roughly one third to one half of the deposition rates measured at the bottom mooring underneath the farm, suggesting that particle deposition beneath the farm consists of more than just oyster biodeposits. While all our deposition rates were within the range of other literature values, ∼10–650 g m⁻² d⁻¹ (Comeau et al., 2014), oyster biodeposition rates can vary widely with season (Mitchell, 2006; Comeau et al., 2014) and oyster size (Haven & Morales-Alamo, 1972), so direct comparisons between farms should be done with that in mind. Deposition at the Away site was particularly high, however in hindsight we suspect this may have been an artifact of the deployment method; the Away mooring was lowered to the bottom from the boat and was repositioned after deployment, which may have created a sediment plume that inflated settlement trap data. We also found that our Away site fell within the deposition footprint of the farm, leaving us unable to parse the relative impact of farm-sourced biodeposits and ambient particle deposition on total particle flux.

Site comparison

The lack of difference in biogeochemical fluxes and macrofaunal community structure between the Farm and Away sites and the relatively low organic matter content indicate that biodeposition is not a major influence at either site, i.e., that regardless of whether the Away site was within the footprint of the farm, deposition was small enough to have minimal impact. Previous studies in other systems have generally found significant differences in biogeochemical fluxes between farm and references sites, however reference sites in these studies tend to be >500 m away from farm sites (Testa et al., 2015; Richard et al., 2007; Giles, Pilditch & Bell, 2006), considerably further than the Away site in this study, which was selected at only ∼90 m from the edge of the farm lease to avoid the flow effects of the constriction in the river downstream. A similar reference site distance was selected by Mallet, Carver & Landry (2006) in a study of an oyster farm in New Brunswick, Canada. There, the researchers also found no difference between sediment geochemistry at reference and farm sites. Additionally, any potential difference in fluxes between our two sites may have been obscured by high within-site variability (Fig. 9). Callier et al. (2009) also found trends, but no significant differences, in biogeochemical fluxes between “unimpacted” reference cores and cores with increasing densities of mussels suspended above them to simulate farm deposition, with high variability within their core treatments.

Like the variability in erosion, microalgal aggregations on the sediment surface may have contributed to variability in fluxes. Microphytobenthos have been shown in several studies to regulate biogeochemical fluxes across the sediment-water interface (Reay, Gallagher & Simmons, 1995; Cerco & Seitzinger, 1997), even when cores are kept in the dark (Sundbäck & Graneli, 1988). Indeed, chl-a measurements from the incubation cores (7.3–22.5 µg g⁻¹ dry sediment) were higher than the highest measured by Watling et al. (2001) (5.24 µg g⁻¹ dry sediment) in sediments from a deeper area of the Damariscotta. Chl-a concentrations were also higher than concentrations of phaeopigments (which indicate degraded chlorophyll often due to grazing), implying abundant living microphytobenthos (e.g., Bianchi, Dawson & Sawangwong, 1988). Chl-a was highly variable between samples, suggesting patchiness that may have contributed to variation in fluxes.

Although there was not a correlation between oxygen flux and DIC (Fig. 9), the ratio of DIC flux to O₂ flux was ∼2, higher than predicted from stoichiometry and observed in high-Arctic sediments (Rysgaard et al., 1998), but similar to other studies in sediments under normoxic conditions and lower than in hypoxic sediments (Lehrter et al., 2012; An & Joye, 2001). This is also in line with a previous study by Newell, Cornwell & Owens (2002) on the influence of biodeposition on sediments. They performed dark incubations of biodeposit-enriched sediment cores with a microalgal mat and found low oxygen influxes relative to DIC efflux, attributing the skewed ratio to the diversion of oxygen into oxidizing reduced compounds (Newell, Cornwell & Owens, 2002).

The “Amended” cores that were enriched with biodeposits prior to incubation showed no differences in nutrient fluxes from unenriched Away site cores, which suggests that these sediments were already saturated with organic matter. The amount of biodeposits added was comparable to the deposition at the farm site over a day (8 mL × 1.4 g mL⁻¹ density of biodeposits spread over the area of a 15 cm diam core is 634 g m⁻²; assuming the biodeposits were loosely packed, this would be an overestimate, so the actual amount added is closer to the 490 g m⁻² d⁻¹ collected in sediment traps under the farm). However, it is important to note that, because of the incubation equipment malfunction, the oxygen flux measurements and nutrient sample collections for the second set of cores did not occur simultaneously, so conclusions about the relationship between oxygen and nutrient fluxes should be drawn with caution.

There was relatively low organic matter content at both sites, supporting the idea that biodeposition is minimal. Surprisingly, sediment at the Away site appeared to be slightly muddier and have higher organic content than the Farm site (Fig. 8). While we expected sediments at the Farm site to be finer and have higher organic content because of increased biodeposition, the potential addition of larger particles from eroding oyster shell to sediment directly beneath the oyster cages could coarsen the sediment and dilute its organic content with large and mostly inorganic shell fragments. However, it is difficult to assess site differences from these data given only one core was taken from each site.

The Away and Farm sites had similar infaunal communities, except for the abundant small nephtyids in some of the Away cores. It is likely that these were juveniles, and the differences in abundance due to patchy recruitment. There is a general absence of recent infaunal community data from this area of the Damariscotta, which makes it difficult to assess from the literature whether the communities at and around the farm are comparable to other locations in the estuary. Additionally, macrofaunal abundance and community assemblage appear very patchy in these sediments and were probably not fully described by our sampling. Also, the cores were stored in a tank with sediment containing large macrofauna from an unrelated experiment, so it is possible that the large nereid and nephtyid worms found when samples were sorted had migrated and established in our cores and were not representative of the natural community at our sites. However, many large nereid jaws were observed in the bulk material while samples were being sorted, as well as many juvenile nephtyids counted, which suggests that the large macrofauna were likely present in the cores when they were collected.

Other than the abundant, small nephtyids at the Away site, the macrofaunal assemblages at both Farm and Away sites mainly consisted of discretely-motile, surface deposit/suspension feeding polychaetes (spionids and flabelligerids) and suspension-feeding infaunal bivalves (myids, mactrids, and pharids). Relative dominance of suspension feeders suggests that flow is slow enough to avoid sediment destabilization but fast enough to maintain sufficient suspended food particle flux (e.g., Widdows et al., 2004). Additionally, surface-deposit feeding may be particularly advantageous where an algal biofilm offers a concentrated high-quality food source at the sediment surface (Decho & Lopez, 1993; Montserrat et al., 2008). Such mats, however, may disadvantage subsurface deposit feeders by impeding sediment-water exchange of oxygen and reduced compounds (e.g., Hansen & Kristensen, 1997). Subsurface deposit feeders contributed less than 20% to Farm and Away total abundance, excluding nephtyid abundance. Though our sampling occurred during a relatively quiescent period of low flow and turbidity, the spring tide several weeks prior (Fig. 7) would have produced high flow rates and sediment erosion potentially stressful for infauna.

Potential flushing of biodeposits

The partial mat failure during the erosion experiments suggests the possibility that these sediments are at least periodically eroded, decreasing the impact of biodeposits on sediments under the farm. The high variability in erosion at the higher shear velocity levels (Fig. 6) suggests sediment erosion occurs suddenly and at shear velocities exceeded only at the highest flow rates recorded by the ADV during the ∼23 h deployment. This is likely because a microphytobenthic mat stabilizes surface sediments, resisting erosion until a point of critical failure; this is consistent with the high chl-a values measured in the incubation cores and with the divers’ observations of patchiness in the mat cover of the sediment surface (Fig. 2). Though only observed in one of our replicate cores, the partial mat failure at the highest tested shear velocity (2.32 cm s⁻¹) falls within the broad range of shear velocities at which biofilms fail. Walker & Grant (2009) performed erosion experiments on sand and sandy mud sediments with patchy cyanobacterial mats beneath mussel long lines, observing mat failure at around 1.5 cm s⁻¹. Grant & Gust (1987), however, examined erosion of purple sulfur bacteria and cyanobacteria biofilms on sand and found critical friction velocity increased from around 1 cm s⁻¹ for clean sands to around 4 cm s⁻¹ on sands with biofilms.

Massive sediment erosion was not observed in our lab experiments, even as the shear velocities tested exceeded the maximum calculated from the measured in situ bottom flow (Fig. 5A). However, the in situ bottom flow was measured only for a short period during neap tide, and higher flow velocities during spring tides would substantially increase shear on the bottom because bed shear stress is proportional to the square of current velocity. Using data from the LOBO station, we predicted that bottom velocity during the time of greatest flow in the spring-neap cycle that summer reached ∼40 cm s⁻¹, driving bottom shear velocities close to 3 cm s⁻¹. This would be well above the shear velocity that resulted in partial mat failure in our erosion experiments. Erosion of surface sediment under the farm during maximum velocities of spring tides probably flushed the deposited sediment away from the farm. The periodic peaks in surface turbidity at spring tides observed in the long-term record are consistent with episodic erosion and flushing of biodeposits. Relatively high eroded mass at two of the three replicates from beneath the farm suggests some flushing of recent biodeposits may have also occurred during the neap tide in which our observations were made (Fig. 6). The poor sorting of sediments at both sites at surface and depth, however, suggests that bulk subsurface sediments do not undergo frequent erosion and deposition events (Fig. 8A, Table S1). Frequent erosion and redeposition tends to result in graded deposition due to gravity sorting and therefore well-sorted sediment profiles (e.g., Morton, 1988). It therefore seems likely that erosion is limited to surface, recently deposited sediments and that these sediments are not redeposited in bulk around the farm.