Jennifer L. Ruesink1*
Katie Houle2
Kristy J. Kroeker3
Brett R. Dumbauld4Fiona C. Boardman1Nathaniel S. Lewis5Brooke A. McIntyre4
Andrew D. Suhrbier2Bobbi Hudson2- 1Department of Biology, University of Washington, Seattle, WA, United States
- 2Pacific Shellfish Institute, Olympia, WA, United States
- 3Department of Ecology and Evolutionary Biology, University of California, Santa Cruz, Santa Cruz, CA, United States
- 4Agricultural Research Service, United States Department of Agriculture, Newport, OR, United States
- 5Coastal Oregon Marine Experiment Station, Oregon State University, Newport, OR, United States
Species persist through resistance and recovery traits induced by natural stressors. Whether human activities exacerbate or moderate the effects of natural stressors is an open question. Because eelgrass (Zostera marina) creates coastal habitat, its response to natural stressors and human activities is of particular management importance. In this study, traits of intertidal eelgrass were examined across cumulative stressors of emersion and oyster culture, including two culture types: oysters grown directly on sediment (ground culture) or supported by gear above the sediment (off-bottom culture). Summer eelgrass was larger above- and below-ground and branched less than in spring, while density and cover were similar seasonally but declined at higher elevation and in ground culture. Eelgrass traits were divided into those related to resistance (larger above-ground size and below-ground storage) and recovery (increased flowering, branching, and rhizome extension). Resistance traits responded additively to intertidal elevation and oyster culture, with above- and below-ground size reduced in all conditions except that rhizome mass was maintained in ground culture. Smaller above-ground size may confer resistance to emersion stress, as a departure from expectations for other stressors. For resilience traits, flowering increased and internode length declined at higher elevation, whereas these traits did not change in ground culture, and off-bottom culture was associated with shorter internodes (additive cumulative stressor) and tended to moderate the enhanced flowering at higher elevations (non-additive, P=0.058). Transitory disturbance in ground culture may reduce eelgrass density with few effects on resilience, whereas off-bottom culture involves longer-term gear placement and trait responses by eelgrass.
1 Introduction
Resilience can occur because a system is resistant to adverse conditions or because it recovers rapidly once favorable conditions return (Nimmo et al., 2015). Across seagrass species, resistance and recovery traits generally show tradeoffs (Kilminster et al., 2015; O’Brien et al., 2018), while they collectively contribute to resilience. Similarly, populations within a species can vary in their ability to resist and recover from environmental stress (Tuya et al., 2021). A wide variety of indicators change when seagrasses are exposed to unfavorable environmental conditions (Roca et al., 2016), which is how we use the term “stressor” here. When stressors exceed resistance and recovery capacities of seagrasses, regime shifts can occur to unvegetated states (Unsworth et al., 2015). For intertidal marine organisms, emersion stress varies as a function of elevation, due to drying out at low tide and/or thermal damage in air (Connell, 1972; Peterson, 1991). Coastal human activities may introduce new stressors or modify existing ones associated with elevation gradients. A key question is the extent to which resilience is maintained under human influence, or alternatively, if anthropogenic activities impair resistance and recovery traits (Cabaço and Santos, 2012; Vieira et al., 2020). This study addressed the response of eelgrass (Zostera marina) to two types of oyster culture practices across an intertidal stress gradient.
How seagrass responds to multiple stressors has been tested for combinations of temperature, salinity, light, and nutrients (Stockbridge et al., 2020; Ostrowski et al., 2021). However, other factors are at play where intertidal seagrass overlaps with aquaculture. In terms of biological interactions, bivalves can compete for space with seagrass but are increasingly recognized for facilitative roles as bivalves filter water, cycle nutrients, and stabilize sediment (Fales et al., 2020; Gagnon et al., 2020). In terms of physical effects, structures in which bivalves are grown can influence light, water flow, and sediment dynamics, and harvest can cause disturbance (Dumbauld et al., 2009; Forrest et al., 2009; Ferriss et al., 2019; Howarth et al., 2022). In the present study region, ground culture refers to the placement of oysters at low density directly on the sediment surface, often mechanically harvested with a boat-towed dredge, while off-bottom culture occurs when oysters are anchored above the sediment surface on lines strung across vertical pipes. Ground culture places oysters at low densities unlikely to compete with eelgrass (Wagner et al., 2012) but periodically removes eelgrass as oysters are harvested (Tallis et al., 2009), acting as a transitory perturbation (Oliver et al., 2015). Off-bottom culture has more gear present but less intense disturbance at harvest, with gear placement generating a rapid-onset perturbation (Oliver et al., 2015). Ferriss et al. (2019) used a meta-analysis to demonstrate reduced density but enhanced growth and reproduction of eelgrass (Zostera spp.) in areas of ground culture, and reduction of density and several traits in off-bottom oyster culture. This result does not address whether bivalve culture effects differ across intertidal stress gradients and whether these two types of aquaculture might differentially affect seagrass traits associated with resilience.
Intertidal zones provide several gradients relevant to seagrass resilience traits. Seagrasses can be light-limited at depth and desiccation-limited at their upper limits, while hydrodynamic energy can act either across or within elevations (Koch, 2001; Yang et al., 2016). These physical factors can require different morphological or demographic strategies at different positions along elevation gradients for seagrass to persist. In terms of size, smaller shoots (=above-ground bundle of leaves) occur with increased water energy and at higher elevations (Keller and Harris, 1966; Boyé et al., 2022). In terms of reproduction, flowering may have a unimodal response to elevation (Olesen et al., 2017). Some responses to environmental gradients, such as larger shoots at depth, may be adaptive for resource (light) acquisition, whereas smaller shoot size may enable closer association with damp sediment and protect from emersion (Shafer et al., 2007). Alternatively, rather than an induction of favorable traits, differential expression of traits across elevational gradients may reflect constraints. Smaller size intertidally could be a result of regular damage (desiccation: (Boese et al., 2005), clipping: (Hernán et al., 2021)), and low light could limit energetic allocation to flowering in subtidal plants (Olesen et al., 2017). Overall, both resistance and recovery traits could involve phenotypic remodeling in ways that generate a morphotype more suited to local conditions, and, particularly for above-ground shoot size, smaller size could improve resistance intertidally. In contrast, across species and in contexts other than intertidal stress, smaller seagrass species are generally considered less stress-resistant than larger species (O’Brien et al., 2018). Because our conclusions are sensitive to this decision to ascribe smaller size as a resistance trait intertidally, we address alternatives in the discussion.
Morphological and reproductive traits that confer resilience (O’Brien et al., 2018) and change with aquaculture (Ferriss et al., 2019) also feature in trait-based frameworks that provide mechanisms for species distributions (Moreira-Saporiti et al., 2023). In part, intraspecific trait variation in Z. marina derives from phenotypic plasticity, which is clearly exhibited in seasonal variation in shoot size and branching (Orth and Moore, 1986; Ruesink et al., 2010). Moreover, some aspects of sexual reproductive investment depend on environmental conditions (Backman and Barilotti, 1976; van Katwijk and van Tussenbroek, 2023). Seasonal expression of eelgrass traits often shifts with latitude and temperature; specifically, there is a delay in peak flowering and biomass at higher latitudes and with cooler water temperatures (Clausen et al., 2014; Blok et al., 2018). Warmer sea surface temperatures have been linked to increased (Potouroglou et al., 2014) and decreased flowering (Thom et al., 2003; Qin et al., 2020). Thus, resilience to small-scale stressors such as emersion, flow, and human activities in the intertidal zone could shift due to broader seasonal or interannual change.
In this study, we addressed a critical gap regarding multiple stressors to seagrass by examining traits of Z. marina across intertidal stress and oyster culture methods, replicated across 12 sites along the US west coast and Salish Sea. We tested how eelgrass traits responded to natural stressors and were modified in ground culture or off-bottom culture of oysters (Figure 1), near the upper elevation limit of eelgrass. Sampling positions spanned a range of eelgrass cover, while tidal elevation was added for potential explanatory power during analysis. We hypothesized:
1. a decline in the amount of eelgrass (density, percent cover) at higher tidal elevations and with oyster culture.
2. trait variation across intertidal gradients providing resistance to emersion (smaller shoot size, increased rhizome mass) and enhanced recovery (increased sexual and asexual reproduction, and rhizome extension).
3. trait variation specific to each oyster culture method, with recovery traits enhanced by transitory disturbance (ground culture) and resistance traits affected by off-bottom culture. These effects could be particularly pronounced at stressful intertidal elevations, generating non-additive statistical effects of elevation and aquaculture (Figure 1).
Figure 1. Response of eelgrass to multiple stressors. (A) Four stressors that can reduce the amount of eelgrass. While density and percent cover may decline due to environmental stress, remodeling could confer resistance (smaller shoot size, larger rhizome mass), and allocation could promote recovery (increased branching, flowering, and expansion). (B, C) Examples show how a trait such as size can decline at higher elevations, with an additional response to aquaculture that is either a main effect (B) or interactive effect (C). (D, E) Effect sizes consistent with multi-stressor responses in (B, C).
2 Materials and methods
2.1 Study design
2.1.1 Environmental conditions at study sites
Eelgrass (Zostera marina) has a broad distribution throughout temperate latitudes in the northern hemisphere (Green and Short, 2003). We studied intertidal eelgrass in a mesotidal region of the NE Pacific Ocean, where the upper limits of the distribution emerge regularly on low tides. At 12 sites in five bays (45.5-48°N; Table 1, Figure 2), eelgrass was sampled across an intertidal gradient, in and outside cultivated beds of oysters (Magallana gigas). Two additional bays in California were sampled in eelgrass meadows outside aquaculture (Figure 2). Mean tidal range spans 1.8 to 2.4 m across sites (tidesandcurrents.noaa.gov), and sampling positions were between -0.6 and +0.6 m relative to mean lower low water (Supplementary Figure S1). Sediment was muddy sand, slightly less muddy at higher elevations (Pearson’s r = -0.28, t(114) = -3.07, P=0.003), corresponding to a typical pattern of grain size increasing at higher elevations (Peterson, 1991). Organic and mud content were strongly correlated (Pearson’s r = 0.81, t(117)=15.08, P<0.0001) (Supplementary Figure S1). Water temperatures were around 10°C in spring at all sites but more variable in summer, ranging from 12.7°C in Tillamook Bay to 19.0°C in Willapa Bay (Supplementary Figure S2).
Table 1. Bays, sites, and sampling positions included in the study design: three levels of aquaculture (none, ground, off-bottom) crossed with three levels of eelgrass (none, sparse, dense).
Figure 2. Map of study sites in (A) Washington, Oregon and California states, USA. Insets of Washington state sites include (B) Samish Bay: PL (Plant), IS (Island); (C) Hood Canal: RP (Rock Point), HH (Hama Hama); (D) Willapa Bay North: CO (Cut Off), WC (West Channel), BC (Bay Center); Willapa Bay South: PO (Port), LI (Long Island), MS (Middle Sands). Oyster grow-out methods are shown in images to the right. In longline culture, clusters of oysters grow on lines strung across stakes <50 cm above the sediment. In flip bag culture, dozens of single oysters grow in mesh bags attached to lines strung across stakes ~70 cm above the sediment, and the bags upend on higher tides due to the floats on the distal side.
Tidal elevation was considered a stressor due to emersion time at each study position. Bathymetry in GIS layers was most fully developed for Willapa Bay, where tidal elevations had been surveyed in 2002 by LiDAR (National Oceanic and Atmospheric Administration) and in 2005 by GPS (US Department of Agriculture). In other bays, publicly-available bathymetry did not align well with personal experience of the duration of low tide. Instead, the time when water flooded at a sampling position was used to infer tidal elevation, based on where the water level was predicted to be at that time at the closest tide station (tidesandcurrents.noaa.gov; Supplementary Table S1). Video cameras were placed at each position, and the time recorded when the water level reached the camera lens. This time was matched to the corresponding time of observed water level at the nearest buoy station. The height of the camera (0.69 m above sediment) was subtracted from this value as a proxy for elevation. At least three observations from spring and summer videos were used to calculate an average elevation for each position.
2.1.2 Oyster aquaculture at study sites
Ground culture was carried out with clusters of oysters at some sites, but single oysters at others. Off-bottom culture occurred on longlines or in flip bags (Table 1). Photographs of the grow-out methods are provided in Figure 2. These culture “treatments” occurred in similar abiotic environments at each site, with three exceptions: MS off-bottom culture at higher elevation, WC organic content higher in ground culture, and BC less organic and mud outside culture (Figure 2, Supplementary Figure S1). Oyster densities on these beds were reported previously and were generally 1-2 year old Pacific oysters at <70 m-2 (Ruesink et al., 2023). One site (MS) had longlines harvested and replanted between the spring and summer sampling, but otherwise beds were sampled in the middle of the crop cycle.
2.1.3 Field data collection at sampling positions
The 12 sites were considered blocked replicates. Sites were selected where three levels of aquaculture occurred within <1 km (none, ground culture, off-bottom culture), and then within each aquaculture level, we specifically selected 3-5 positions with different amounts of eelgrass, ranging from none to the highest cover available (Table 1). A “position” required at least 1000 m2 with a similar cover of eelgrass, within a given grow-out method or with no aquaculture present. Although habitat and eelgrass traits were measured in multiple quadrats, position is the true replicate. Where possible, eelgrass was measured along a 100-m transect, but for smaller areas, in two 50-m transects. Along each 100-m length, 20 quadrats (0.25 m2) were spaced at 5-m intervals. Vegetative and flowering shoots were counted in each quadrat. Also, at Hood Canal sites, and for longline culture at PL where no eelgrass was collected, maximum length of five terminal shoots was measured per quadrat. At all other sites, up to 5 terminal shoots were collected with attached rhizome from each of 10 quadrats. Each position was sampled in spring (Apr-Jun) and summer (Jul-Aug). Two seasons were specifically included in the design because of typical changes in traits and density between the beginning and end of the growing season (Clausen et al., 2014; Blok et al., 2018; Ruesink et al., 2022), which also occurred in eelgrass meadows in the study region (Supplementary Figure S3).
2.1.4 Eelgrass trait measurements
Collected shoots were measured for maximum length, sheath length, sheath width, and the length of the most recent fully-extended internode on the rhizome (longest internode length within two nodes; Supplementary Figure S4). Above-ground material (leaves and sheaths) was dried (60°C for at least 5 days) and weighed. A section of rhizome up to 7 cm in length was dried and weighed separately, then standardized to mass per unit length. Any branches that had formed within the sheath or at the two most recent nodes on the rhizome were counted (and then removed before determining biomass of the terminal shoot).
For cover and density data, 4738 quadrats were measured in two seasons across 12 sites, whereas shoot length and flowering data came only from quadrats in positions selected to include eelgrass (3178 quadrats). Other traits were summarized seasonally across all terminal, non-flowering shoots collected from each position across 10 sites (4868 shoots). Prior to statistical analysis, quadrats within each position were averaged for cover, density, flowering frequency, and shoot length; shoots collected from each position were averaged for other traits. Based on true replicates (positions), sample size for each analysis was 118 for cover and density across 12 sites, 75 for flowering and maximum length across 12 sites, and 63 for traits of shoots collected across 10 sites.
2.2 Data analysis
2.2.1 Amount of eelgrass across stressors
At each site, positions were deliberately selected where no eelgrass was present and in sparse and dense eelgrass (Table 1), so the analysis addressed whether small differences in tidal elevation could underlie some of this variability. Also, the design would not be expected to detect differences in the amount of eelgrass across aquaculture “treatments” unless the upper range of densities was altered on aquaculture beds. Nevertheless, we analyzed whether the amount of eelgrass (percent cover in quadrats at low tide, shoot density) responded to the following factors: season, elevation, and oyster culture (none, ground, off-bottom). Initially we included all two- and three-way interactions, but, since statistically-significant interactions of season with other factors did not occur, and we were primarily interested in cumulative stressors, generalized linear mixed models (GLMM) were built with season, elevation, and culture as main effects, plus an elevation x culture interaction. (Zuur et al., 2009) recommend a model comparison approach to determine the appropriate random effect to include (i.e., intercept, intercept and slope, no random effect). Generally, models fit with a random slope failed to converge, and we elected to use site as a random intercept in all models to appropriately account for the study design, which was blocked by site. This random site effect accounts for sites studied in different years as well. Residuals were examined to select an appropriate data distribution (glmmTMB (Brooks et al., 2017), DHARMa packages (Hartig, 2022)) in R (R Core Team, 2023). Datasets of density and cover, which contained many zeros, were fit assuming a Tweedie distribution (log-link). Tweedie falls within the general class of exponential dispersion models and was preferred over zero-inflated models that require interpretation of two outputs: the probability of presence, and the amount when present. Final model results are displayed as effect sizes and 95% confidence intervals for main effects of season, elevation, oyster culture (ground and off-bottom compared to no culture), and the culture x elevation interaction (sjPlot package (Lüdecke, 2023)). Explanatory power of each GLMM was calculated as McFadden’s pseudo-R2 (Shtatland et al., 2002).
2.2.2 Multivariate trait analysis
Exploratory analyses for eelgrass trait correlations were performed by principal components analysis, normalizing traits by their mean and standard deviation. Specifically, these traits were maximum length, sheath length, sheath width, internode length, above-ground biomass per shoot, below-ground biomass per cm of rhizome, recent branching, and flowering frequency. Each position was a sample, visualized by season and site. Density and cover were not included in this multivariate trait analysis, as these are attributes of the position, not individual-level traits, and additionally positions were deliberately selected to span a large range of amount of eelgrass at each site and in each level of culture.
2.2.3 Traits of eelgrass across stressors
In a multiple-stressors framework, we analyzed two resistance traits (above-ground maximum length, below-ground biomass per cm of rhizome) and three recovery traits (flowering frequency, recent branching, internode length). Because shoots were not collected from two of 12 sites, branching, internode length, and rhizome mass were restricted to samples from 10 sites. GLMMs were generated for each response variable, with season, elevation, and culture as main effects, the elevation x culture interaction, and site as a random effect. Data distributions for morphological variables were assumed Gaussian, and for flowering and branching with many zeros were tweedie (log-link), following inspection of residuals (Supplementary Table S2). For completeness, we also set up models with season, level of eelgrass (sparse, dense), culture, and the eelgrass x culture interaction as predictor variables for eelgrass traits, which aligned with the a priori study design that did not account for tidal elevation. Eelgrass density (sparse, dense) always had poorer relationships with traits than did tidal elevation (Supplementary Table S2).
3 Results
3.1 Amount of eelgrass
Density and cover of eelgrass declined at higher tidal elevations, consistent with intertidal stress reducing the amount of eelgrass. Effect sizes from oyster culture also tended to be negative, significantly so in ground culture (Figure 3). Elevation x culture interactions were not significant: cumulative stressors were additive (Supplementary Table S2). In analyses of both amount and traits of Z. marina, McFadden’s pseudo-R2 values were 0.11-0.31 (corresponding to fixed effects), except lower R2 values for density and below-ground biomass where site differences were pronounced (site-level variation was 0.30 and 0.31 of total variation in density and below-ground biomass).
Figure 3. Percent cover of eelgrass (Zostera marina) in (A) spring and (B) summer across intertidal elevation and oyster aquaculture (Ground and Off-bottom). Lines and 95% confidence intervals are simple linear regressions across sample positions in each aquaculture type, which does not account for site effects incorporated in statistical models. Effect sizes for (C) percent cover and (D) shoot density from three-factor mixed effect model results. The neutral effect line is placed at 1 due to log-link in error structure. Positive effects are right of the neutral line, negative effects to left. Error bars show 95% confidence intervals of effect sizes; error bars that cross the neutral effect line indicate non-significant effects at α=0.05. Season denotes summer relative to spring, Ground and Off-bottom are in relation to outside aquaculture. N= 236 (118 positions across 12 sites, each measured in spring and summer).
3.2 Multivariate trait analysis
The first two principal components explained 68% of total variation in eight eelgrass traits across positions and seasons (Figure 4). Above-ground size variables loaded strongly on PC1. Below-ground variables, specifically internode length and rhizome mass, were less strongly aligned with this size axis. Flowering frequency and recent branching loaded in opposite directions on PC2. Spring samples tended to fall in the lower right-hand quadrant of the component scores, in keeping with smaller shoots that branch frequently in that season. Summer shoots were generally larger with greater flowering frequencies (Figure 4). This multivariate analysis confirmed the importance of sampling in two seasons, while motivating a choice of five traits for univariate analysis: maximum length and below-ground biomass as resistance traits, and flowering, branching, and internode length as recovery traits.
Figure 4. Principal components analysis of eight traits of eelgrass (Zostera marina) at 10 sites in Washington and Oregon, USA. Site codes are in Table 1. PC1 explained 52% and PC2 explained 16% of total variation. Points are sampling positions (60 in spring, 63 in summer). SHL, sheath length; MAX, maximum shoot length; SHW, sheath width; AGB, above-ground shoot biomass; BGB, rhizome biomass as dry mass per length; INT, internode length; FLO, flowering frequency; BRA, recent branches.
3.3 Resistance traits
Shoot length declined at higher elevations and with both oyster culture methods (Figure 5A). Emersion and culture stressors together were estimated to reduce shoot length by one-third to one-half, based on model estimates (Supplementary Table S2). Biomass per rhizome length was reduced at higher elevation and in off-bottom culture, but not in ground culture (Figure 5B). No significant interaction effects were found between the two stressors. Eelgrass shoots in summer were larger both above- and below-ground than in spring (Supplementary Table S2).
Figure 5. Effect sizes for resistance traits of eelgrass (Zostera marina) across intertidal elevation and oyster aquaculture (Ground and Off-bottom). (A) Maximum shoot length (bundle of leaves) of terminal vegetative shoots at 12 sites (75 positions, 2 seasons). (B) Dry biomass of a standardized length of rhizome for shoots collected from 10 sites (63 positions, 2 seasons). The neutral line is placed at 0 due to Gaussian error structure. Display of effect sizes follows Figure 3.
3.4 Recovery traits
Flowering frequency increased at higher elevations, but the increase at higher elevation tended to be moderated in off-bottom culture (negative interaction effect, P=0.058; Supplementary Table S2). Flowering had a positive effect size for ground culture, although not significant (Figure 6A). Branching frequency averaged 0.16/node in spring and 0.08/node in summer, and season was the only significant predictor (Figure 6B, Supplementary Table S2). Internode length was shorter at higher elevations and in off-bottom culture (Figure 6C, Supplementary Table S2).
Figure 6. Effect sizes for recovery traits of eelgrass (Zostera marina) across intertidal elevation and oyster aquaculture (Ground and Off-bottom). (A) Flowering frequency based on quadrats at 12 sites (75 positions, 2 seasons). (B) Branching and (C) internode length of shoots collected from 10 sites (63 positions in 2 seasons). The neutral line is placed at 0 for Gaussian error structure and at 1 for flowering and branching with log-link error structure. Display of effect sizes follows Figure 3.
4 Discussion
Numerical and cover reductions reliably indicate stress for seagrass (Roca et al., 2016). In this study, the amount of eelgrass consistently declined at higher elevations, and oyster culture added stress, even though we deliberately selected positions in and out of culture across a full range of eelgrass densities. Intraspecific trait variation revealed a more complex set of responses to stressors, with some cases showing eroded resilience while others maintained resistance or recovery traits in oyster aquaculture. Eelgrass in ground culture was smaller above- but not below-ground. In off-bottom culture, eelgrass size decreased both above- and below-ground, and while branching did not change, recovery traits of flowering and rhizome extension diminished in some parts of the tidal elevation gradient.
Our field study design, where stressor treatments occurred at each site (blocked design), highlighted environmentally-modulated trait changes. However, due to low replication of treatments within sites, we could not explore how tidal elevation responses might differ under varying culture methods at each site (including e.g., singles vs. clusters, or longlines vs. flip bags). Overall, the impact of aquaculture on resilience traits of eelgrass depended on the oyster grow-out method (Figures 1, 5, 6). This general result aligns with a meta-analysis showing that both ground culture and off-bottom culture reduce the amount of eelgrass, with only off-bottom culture affecting reproduction (Ferriss et al., 2019). Multiple stressors often act additively on seagrasses (Stockbridge et al., 2020), and this observational study was no exception. Statistically, only flowering showed an elevation x culture interaction, with large but not significant effect size (P=0.058), in which the enhanced flowering typical of higher elevations was induced less in off-bottom culture (Figure 6A).
Multivariate trait analyses in seagrasses commonly visualize correlated traits (e.g (Nomme and Harrison, 1991)). When density and traits are both included, smaller shoots typically occur at higher density, with a second axis corresponding to the amount of eelgrass (i.e., cover, mass per area, or leaf area index; (Olesen and Sand Jensen, 1994; Vieira et al., 2018; Boyé et al., 2022; Duffy et al., 2022)). In our analysis (Figure 4), we focused on individual traits rather than population-level characteristics (density, cover). Both multivariate and individual trait analyses showed a seasonal shift from small, rapidly-branching shoots in spring to larger shoots in summer, when flowering became more evident (Figure 4, 5, 6, Supplementary Figure S3). Our results confirm that longer shoots share covarying traits like wider and longer sheaths, faster leaf and rhizome extension, and larger rhizomes (Ruesink et al., 2018; Boyé et al., 2022). This covariation makes it notable that shoot length but not rhizome mass differed in ground culture (Figure 5), a functionally-important decoupling of above- and below-ground responses. Smaller size of shoots above-ground may protect from further damage, but smaller storage below ground could undermine eelgrass’s ability to survive stressful conditions. Leaf damage may drive phenotypic plasticity in above-ground size, as subsequent leaves emerge smaller and less exposed to desiccation (Ruesink et al., 2012).
Tidal elevation strongly influenced both the amount and traits of eelgrass. Cover and density declined at higher elevations, consistent with sampling near the upper limit of eelgrass (Figure 3). Both above- and below-ground size declined at higher elevation, generally reducing resistance traits. Recovery traits were either induced (flowering) or reduced (internode length) at higher elevation, and branching also had a negative mean response to elevation. Although flowering is a resilience trait, recovery depends on successful seed set and germination, which we did not investigate. Shoot size has been shown repeatedly to decline at higher intertidal elevations (Keller and Harris, 1966; Boese and Robbins, 2008; Ruesink et al., 2012). Increased flowering at higher intertidal elevations has been reported in Willapa Bay and other coastal estuaries (Boese and Robbins, 2008; Ruesink et al., 2012). These responses, while common intertidally, are not universal. Flowering frequency may be unresponsive to depth (Yang et al., 2013; von Staats et al., 2021; Lekammudiyanse et al., 2024). Subtidal eelgrass studies report varying trait patterns with depth, such as taller shoots but unimodal relationships for amount of eelgrass and investment in flowering (Olesen et al., 2017). In Z. marina, long-term shading can reduce but sediment nutrients may increase flowering frequency (Backman and Barilotti, 1976; Johnson et al., 2017). Since our study occurred in the low-intertidal zone, we did not expect light limitations, unlike in subtidal eelgrass or under floating oyster bags where growth and flowering decline (Dennison and Alberte, 1985; Skinner et al., 2013), although wide gear spacing allows persistent seagrass (Ferretto et al., 2022). The timescale of our study cannot address interannual flowering variation, which can be substantial. We found <1% flowering in Hood Canal and Samish Bay in this study, but these bays have previously been reported to have up to 7% and 10% flowering, respectively (Yang et al., 2013).
While elevation had strong statistical effects, emersion is not definitively the most important abiotic stressor for eelgrass, especially since we did not stratify samples across tidal elevations. Sediment properties, which reflect the overlying hydrodynamic regime, could also influence eelgrass traits through nutrient availability, phytotoxins, or anchoring (Koch, 2001). Lower-elevation sediments were generally muddier, which could drive trait responses if nutrient-rich muddy sediments promote growth or sandy sediments indicative of higher energy induce smaller shoots with low drag (Baden and Boström, 2001; Yang et al., 2013). Internode length might respond to local resource conditions as a foraging strategy to remain in high-resource patches (short internodes) or access new resources (long internodes) (De Kroon and Hutchings, 1995). We did not capture a range of sediment conditions that appeared stressful, so it was also not possible to address whether oyster culture could moderate this stressor.
Oyster aquaculture along the US west coast occurs on broad tidal flats overlapping with eelgrass at low intertidal elevations (Dumbauld and McCoy, 2015) where eelgrass already experiences strong environmental gradients of stress and seasonal light and thermal constraints. Our findings suggest that oyster aquaculture acts additively with natural stressors in most cases, with non-additive (marginal significance) effects observed in higher-elevation off-bottom culture where flowering decreased (Figure 1). Resistance and recovery, especially if induced by stress, have been essential for eelgrass to persist in the variable intertidal environment. Declining resilience traits, particularly when eelgrass amounts are low, raise concerns. Future work can therefore address not only aquaculture practices that maintain co-located seagrass (e.g (Ferretto et al., 2022)), but also how they may enhance eelgrass resistance and recovery traits.
Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: Mendeley Data DOI: 10.17632/d8vgk47fnz; https://data.mendeley.com/datasets/d8vgk47fnz.
Author contributions
JR: Writing – review & editing, Writing – original draft, Visualization, Supervision, Methodology, Investigation, Funding acquisition, Formal analysis, Data curation, Conceptualization. KH: Writing – review & editing, Investigation, Data curation. KK: Writing – review & editing, Validation, Supervision, Investigation. BD: Writing – review & editing, Resources, Investigation, Data curation. FB: Writing – review & editing, Investigation, Data curation. NL: Writing – review & editing, Investigation. BM: Writing – review & editing, Investigation. AS: Writing – review & editing, Investigation. BH: Writing – review & editing, Project administration, Funding acquisition, Conceptualization.
Funding
The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This publication was prepared under Pacific Shellfish Institute grant (#20-31G) with the Pacific States Marine Fisheries Commission (Award #NA18NMF4720007). The statements, findings, conclusions, and recommendations are those of the authors and do not necessarily reflect the views of the National Oceanic and Atmospheric Administration or the US Government. Collaborative work was partially funded by USDA-ARS (CRIS project 2072-63000-005-00D).
Acknowledgments
Field work was carried out with help from WW Hull, E Subbotin, AC Trimble, Natalie Sahli, Kalloway Page, Evie Fagergren, Aimee Christy, Isabel Platten, S Lummis, B Walker, S Myers, E Herman, Joseph Brockman, Julianne Merrill, and Austin Williams. Aspen Katla helped compile water temperature records, and Kalloway Page extracted tidal elevations from water levels. We appreciate site access granted by Taylor Shellfish, Pacific Shellfish, Northern Oyster Company, Jolly Roger Oysters, Goose Point Oysters, R&B Oyster Company, Hama Oyster Company, Rock Point Oyster, and Hog Island Oyster Company. AI was applied to the discussion section as part of the editing process to reduce use of passive voice.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmars.2024.1427595/full#supplementary-material
References
Backman T. W., Barilotti D. C. (1976). Irradiance reduction: effects on standing crops of the eelgrass Zostera marina in a coastal lagoon. Mar. Biol. 34, 33–44. doi: 10.1007/BF00390785
Baden S. P., Boström C. (2001). “The leaf canopy of seagrass beds: faunal community structure and function in a salinity gradient along the Swedish coast,” in Ecological comparisons of sedimentary shores. Ecological Studies, vol. 151 . Ed. Reise K. (Berlin, Heidelberg: Springer-Verlag), 213–236.
Blok S. E., Olesen B., Krause-Jensen D. (2018). Life history events of eelgrass Zostera marina populations across gradients of latitude and temperature. Mar. Ecol. Prog. Ser. 590, 79–93. doi: 10.3354/meps12479
Boese B. L., Robbins B. D. (2008). Effects of erosion and macroalgae on intertidal eelgrass (Zostera marina) in a northeastern Pacific estuary (USA). Botanica Marina 51, 247–257. doi: 10.1515/BOT.2008.034
Boese B. L., Robbins B. D., Thursby G. (2005). Desiccation is a limiting factor for eelgrass (Zostera marina L.) distribution in the intertidal zone of a northeastern Pacific (USA) estuary. Botanica Marina 48, 274–283. doi: 10.1515/BOT.2005.037
Boyé A., Gauthier O., Becheler R., Le Garrec V., Hily C., Maguer M., et al. (2022). Drivers and limits of phenotypic responses in vulnerable seagrass populations: Zostera marina in the intertidal. J. Ecol. 110, 144–161. doi: 10.1111/1365-2745.13791
Brooks M. E., Kristensen K., Benthem K. J., van, Magnusson A., Berg C. W., Nielsen A., et al. (2017). glmmTMB balances speed and flexibility among packages for zero-inflated generalized linear mixed modeling. R J. 9, 378–400. doi: 10.32614/RJ-2017-066
Cabaço S., Santos R. (2012). Seagrass reproductive effort as an ecological indicator of disturbance. Ecol. Indic. 23, 116–122. doi: 10.1016/j.ecolind.2012.03.022
Clausen K. K., Krause-Jensen D., Olesen B. (2014). Seasonality of eelgrass biomass across gradients in temperature and latitude. Mar. Ecol. Prog. Ser. 506, 71–86. doi: 10.3354/meps10800
Connell J. H. (1972). Community interactions on marine rocky intertidal shores. Annu. Rev. Ecol. System. 3, 169–192. doi: 10.1146/annurev.es.03.110172.001125
De Kroon H., Hutchings M. J. (1995). Morphological plasticity in clonal plants: the foraging concept reconsidered. J. Ecol. 83, 143–152. doi: 10.2307/2261158
Dennison W. C., Alberte R. S. (1985). Role of daily light period in the depth distribution of Zostera marina (eelgrass). Mar. Ecol. Prog. Ser. 25, 51–61. doi: 10.3354/meps025051
Duffy J. E., Stachowicz J. J., Reynolds P. L., Hovel K. A., Jahnke M., Sotka E. E., et al. (2022). A Pleistocene legacy structures variation in modern seagrass ecosystems. Proc. Natl. Acad. Sci. 119 (32), e2121425119. doi: 10.1073/pnas.2121425119
Dumbauld B. R., McCoy L. M. (2015). Effect of oyster aquaculture on seagrass Zostera marina at the estuarine landscape scale in Willapa Bay, Washington (USA). Aquacult. Environ. Interact. 7, 29–47. doi: 10.3354/aei00131
Dumbauld B. R., Ruesink J. L., Rumrill S. S. (2009). The ecological role of bivalve shellfish aquaculture in the estuarine environment: A review with application to oyster and clam culture in West Coast (USA) estuaries. Aquaculture 290, 196–223. doi: 10.1016/j.aquaculture.2009.02.033
Fales R. J., Boardman F. C., Ruesink J. L. (2020). Reciprocal interactions between bivalve molluscs and seagrass: A review and meta-analysis. J. Shellf. Res. 39 (3), 547–562. doi: 10.2983/035.039.0305
Ferretto G., Vergés A., Poore A. G. B., Gribben P. E., Glasby T. M. (2022). Floating bags have the potential to minimise oyster farming impacts on Posidonia australis seagrass meadows. Aquaculture 560, 738594. doi: 10.1016/j.aquaculture.2022.738594
Ferriss B. E., Conway-Cranos L. L., Sanderson B. L., Hoberecht L. (2019). Bivalve aquaculture and eelgrass: A global meta-analysis. Aquaculture 498, 254–262. doi: 10.1016/j.aquaculture.2018.08.046
Forrest B. M., Keeley N. B., Hopkins G. A., Webb S. C., Clement D. M. (2009). Bivalve aquaculture in estuaries: Review and synthesis of oyster cultivation effects. Aquaculture 298, 1–15. doi: 10.1016/j.aquaculture.2009.09.032
Gagnon K., Rinde E., Bengil E. G. T., Carugati L., Christianen M. J. A., Danovaro R., et al. (2020). Facilitating foundation species: The potential for plant–bivalve interactions to improve habitat restoration success. J. Appl. Ecol. 57, 1161–1179. doi: 10.1111/1365-2664.13605
Green E. P., Short F. T. (2003). World Atlas of Seagrasses. Eds. Green E. P., Short F. T. (Berkeley USA: University of California Press).
Hartig F. (2022). DHARMa: Residual Diagnostics for Hierarchical (Multi-Level/Mixed) Regression Models. Available online at: https://CRAN.R-project.org/package=DHARMa (Accessed April 25, 2023).
Hernán G., Ortega M. J., Henderson J., Alós J., Boyer K., Cimon S., et al. (2021). Latitudinal variation in plant defence against herbivory in a marine foundation species does not follow a linear pattern: The importance of resource availability. Global Ecol. Biogeogr. 30, 220–234. doi: 10.1111/geb.13217
Howarth L., Lewis-McCrea L., Kellogg L., Apostolaki E., Reid G. (2022). Aquaculture and eelgrass Zostera marina interactions in temperate ecosystems. Aquacult. Environ. Interact. 14, 15–34. doi: 10.3354/aei00426
Johnson A. J., Moore K. A., Orth R. J. (2017). The influence of resource availability on flowering intensity in Zostera marina (L.). J. Exp. Mar. Biol. Ecol. 490, 13–22. doi: 10.1016/j.jembe.2017.02.002
Keller M., Harris S. W. (1966). The growth of eelgrass in relation to tidal depth. J. Wildl. Manage. 30, 280–285. doi: 10.2307/3797815
Kilminster K., McMahon K., Waycott M., Kendrick G. A., Scanes P., McKenzie L., et al. (2015). Unravelling complexity in seagrass systems for management: Australia as a microcosm. Sci. Total Environ. 534, 97–109. doi: 10.1016/j.scitotenv.2015.04.061
Koch E. W. (2001). Beyond light: physical, geological, and geochemical parameters as possible submersed aquatic vegetation habitat requirements. Estuaries 24, 1–17. doi: 10.2307/1352808
Lekammudiyanse M. U., Saunders M. I., Flint N., Irving A., Aiken C., Clark D. E., et al. (2024). Environmental drivers of flowering in the genus Zostera and spatio-temporal variability of Zostera muelleri flowering in Australasia. Aquat. Conserv.: Mar. Freshw. Ecosyst. 34, e4068. doi: 10.1002/aqc.4068
Lüdecke D. (2023). sjPlot: Data Visualization for Statistics in Social Science. Available online at: https://CRAN.R-project.org/package=sjPlot (Accessed April 25, 2023).
Moreira-Saporiti A., Teichberg M., Garnier E., Cornelissen J. H. C., Alcoverro T., Björk M., et al. (2023). A trait-based framework for seagrass ecology: Trends and prospects. Front. Plant Sci. 14, 1088643. doi: 10.3389/fpls.2023.1088643
Nimmo D. G., Mac Nally R., Cunningham S. C., Haslem A., Bennett A. F. (2015). Vive la résistance: reviving resistance for 21st century conservation. Trends Ecol. Evol. 30, 516–523. doi: 10.1016/j.tree.2015.07.008
Nomme K. M., Harrison P. G. (1991). A multivariate comparison of the seagrasses Zostera marina and Zostera japonica in monospecific versus mixed populations. Can. J. Bot. 69, 1984–1990. doi: 10.1139/b91-249
O’Brien K. R., Waycott M., Maxwell P., Kendrick G. A., Udy J. W., Ferguson A. J. P., et al. (2018). Seagrass ecosystem trajectory depends on the relative timescales of resistance, recovery and disturbance. Mar. pollut. Bull. 134, 166–176. doi: 10.1016/j.marpolbul.2017.09.006
Olesen B., Krause-Jensen D., Christensen P. B. (2017). Depth-related changes in reproductive strategy of a cold-temperate Zostera marina meadow. Estuaries Coasts 40, 553–563. doi: 10.1007/s12237-016-0155-4
Olesen B., Sand Jensen K. (1994). Biomass-density patterns in the temperate seagrass Zostera marina. Mar. Ecol. Prog. Ser. 109, 283–291. doi: 10.3354/meps109283
Oliver T. H., Heard M. S., Isaac N. J. B., Roy D. B., Procter D., Eigenbrod F., et al. (2015). Biodiversity and resilience of ecosystem functions. Trends Ecol. Evol. 30, 673–684. doi: 10.1016/j.tree.2015.08.009
Orth R. J., Moore K. A. (1986). Seasonal and year-to-year variations in the growth of Zostera marina L. (eelgrass) in the lower Chesapeake Bay. Aquat. Bot. 24, 335–341. doi: 10.1016/0304-3770(86)90100-2
Ostrowski A., Connolly R. M., Sievers M. (2021). Evaluating multiple stressor research in coastal wetlands: A systematic review. Mar. Environ. Res. 164, 105239. doi: 10.1016/j.marenvres.2020.105239
Peterson C. H. (1991). Intertidal zonation of marine invertebrates in sand and mud. Am. Sci. 79, 236–249. doi: https://www.jstor.org/stable/29774371
Potouroglou M., Kenyon E. J., Gall A., Cook K. J., Bull J. C. (2014). The roles of flowering, overwinter survival and sea surface temperature in the long-term population dynamics of Zostera marina around the Isles of Scilly, UK. Mar. pollut. Bull. 83, 500–507. doi: 10.1016/j.marpolbul.2014.03.035
Qin L.-Z., Kim S. H., Song H.-J., Kim H. G., Suonan Z., Kwon O., et al. (2020). Long-term variability in the flowering phenology and intensity of the temperate seagrass Zostera marina in response to regional sea warming. Ecol. Indic. 119, 106821. doi: 10.1016/j.ecolind.2020.106821
R Core Team (2023). R: A Language and Environment for Statistical Computing (Vienna, Austria: R Foundation for Statistical Computing). Available at: https://www.r-project.org/.
Roca G., Alcoverro T., Krause-Jensen D., Balsby T. J. S., van Katwijk M. M., Marbà N., et al. (2016). Response of seagrass indicators to shifts in environmental stressors: A global review and management synthesis. Ecol. Indic. 63, 310–323. doi: 10.1016/j.ecolind.2015.12.007
Ruesink J., Briones Ortiz B., Mawson C., Boardman F. (2022). Tradeoffs in life history investment of eelgrass Zostera marina across estuarine intertidal conditions. Mar. Ecol. Prog. Ser. 686, 61–70. doi: 10.3354/meps14000
Ruesink J. L., Fitzpatrick J. P., Dumbauld B. R., Hacker S. D., Trimble A. C., Wagner E. L., et al. (2012). Life history and morphological shifts in an intertidal seagrass following multiple disturbances. J. Exp. Mar. Biol. Ecol. 424–425, 25–31. doi: 10.1016/j.jembe.2012.05.002
Ruesink J. L., Hong J.-S., Wisehart L., Hacker S. D., Dumbauld B. R., Hessing-Lewis M., et al. (2010). Congener comparison of native (Zostera marina) and introduced (Z. japonica) eelgrass at multiple scales within a Pacific Northwest estuary. Biol. Invasions 12 (6), 1773–1790. doi: 10.1007/s10530-009-9588-z
Ruesink J. L., Houle K., Beck E., Boardman F. C., Suhrbier A., Hudson B. (2023). Intertidal grow-out technique, not eelgrass (Zostera marina), influences performance of Pacific oysters (Magallana gigas). Aquacult. Res. 2023, 1–13. doi: 10.1155/2023/6621043
Ruesink J. L., Stachowicz J. J., Reynolds P. L., Boström C., Cusson M., Douglass J., et al. (2018). Form-function relationships in a marine foundation species depend on scale: A shoot to global perspective from a distributed ecological experiment. Oikos 127, 364–374. doi: 10.1111/oik.04270
Shafer D. J., Sherman T. D., Wyllie-Echeverria S. (2007). Do desiccation tolerances control the vertical distribution of intertidal seagrasses? Aquat. Bot. 87, 161–166. doi: 10.1016/j.aquabot.2007.04.003
Shtatland E. S., Kleinman K., Cain E. M. (2002). One more time about R2 measures of fit in logistic regression - Northeast SAS User’s Group (NESUG 15). Available online at: https://www.lexjansen.com/nesug/nesug02/st/st004.pdf (Accessed October 16, 2024).
Skinner M., Courtenay S., McKindsey C. (2013). Reductions in distribution, photosynthesis, and productivity of eelgrass Zostera marina associated with oyster Crassostrea virginica aquaculture. Mar. Ecol. Prog. Ser. 486, 105–119. doi: 10.3354/meps10345
Stockbridge J., Jones A. R., Gillanders B. M. (2020). A meta-analysis of multiple stressors on seagrasses in the context of marine spatial cumulative impacts assessment. Sci. Rep. 10, 11934. doi: 10.1038/s41598-020-68801-w
Tallis H. M., Ruesink J. L., Dumbauld B., Hacker S., Wisehart L. M. (2009). Oysters and aquaculture practices affect eelgrass density and productivity in a Pacific Northwest estuary. J. Shellf. Res. 28, 251–261. doi: 10.2983/035.028.0207
Thom R. M., Borde A. B., Rumrill S., Woodruff D. L., Williams G. D., Southard J. A., et al. (2003). Factors influencing spatial and annual variability in eelgrass (Zostera marina L.) meadows in Willapa Bay, Washington, and Coos Bay, Oregon, estuaries. Estuaries 26, 1117–1129. doi: 10.1007/BF02803368
Tuya F., Fernández-Torquemada Y., del Pilar-Ruso Y., Espino F., Manent P., Curbelo L., et al. (2021). Partitioning resilience of a marine foundation species into resistance and recovery trajectories. Oecologia 196, 515–527. doi: 10.1007/s00442-021-04945-4
Unsworth R. K. F., Collier C. J., Waycott M., Mckenzie L. J., Cullen-Unsworth L. C. (2015). A framework for the resilience of seagrass ecosystems. Mar. pollut. Bull. 100, 34–46. doi: 10.1016/j.marpolbul.2015.08.016
van Katwijk M. M., van Tussenbroek B. I. (2023). Facultative annual life cycles in seagrasses. Plants 12, 2002. doi: 10.3390/plants12102002
Vieira V. M. N. C. S., Lopes I. E., Creed J. C. (2018). The biomass–density relationship in seagrasses and its use as an ecological indicator. BMC Ecol. 18, 44. doi: 10.1186/s12898-018-0200-1
Vieira R., Martin A., Engelen A. H., Thomsen M. S., Arenas F. (2020). Interactive effects of co-occurring anthropogenic stressors on the seagrass, Zostera noltei. Ecol. Indic. 109, 105780. doi: 10.1016/j.ecolind.2019.105780
von Staats D. A., Hanley T. C., Hays C. G., Madden S. R., Sotka E. E., Hughes A. R. (2021). Intra-meadow variation in seagrass flowering phenology across depths. Estuaries Coasts 44, 325–338. doi: 10.1007/s12237-020-00814-0
Wagner E., Dumbauld B. R., Hacker S. D., Trimble A. C., Wisehart L. M., Ruesink J. L. (2012). Density-dependent effects of an introduced oyster, Crassostrea gigas, on a native intertidal seagrass, Zostera marina. Mar. Ecol. Prog. Ser. 468, 149–160. doi: 10.3354/meps09952
Yang S., HilleRisLambers J., Ruesink J. L. (2016). Reversal of intraspecific interactions by an ecosystem engineer leads to variable seedling success along a stress gradient. Mar. Ecol. Prog. Ser. 543, 163–171. doi: 10.3354/meps11561
Yang S., Wheat E. E., Horwith M. J., Ruesink J. L. (2013). Relative impacts of natural stressors on life history traits underlying resilience of intertidal eelgrass (Zostera marina L.). Estuaries Coasts 36, 1006–1013. doi: 10.1007/s12237-013-9609-0
Keywords: Zostera marina, stressors, disturbance, perturbation, tidal elevation, eelgrass, resistance, recovery
Citation: Ruesink JL, Houle K, Kroeker KJ, Dumbauld BR, Boardman FC, Lewis NS, McIntyre BA, Suhrbier AD and Hudson B (2024) Intraspecific variation in resilience traits of eelgrass across intertidal stress gradients and oyster aquaculture methods. Front. Mar. Sci. 11:1427595. doi: 10.3389/fmars.2024.1427595
Received: 04 May 2024; Accepted: 07 October 2024;
Published: 01 November 2024.
Edited by:
Rochelle Diane Seitz, College of William & Mary, United StatesReviewed by:
Aurélien Boyé, Institut Français de Recherche pour l’Exploitation de la Mer (IFREMER), FranceFederica Montesanto, University of Florida, United States
Copyright © 2024 Ruesink, Houle, Kroeker, Dumbauld, Boardman, Lewis, McIntyre, Suhrbier and Hudson. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Jennifer L. Ruesink, ruesink@uw.edu