The role of short-term temperature variability and light in shaping the phenology and characteristics of seagrass beds

Abstract

Seagrass beds inhabit highly heterogeneous temperature regimes that characterize the marine nearshore. Temperature directly influences seagrasses and also provides indirect information on other ecologically relevant environmental variables. Multiple temperature processes operate on seasonal and sub-seasonal timescales (i.e., hours to months) and include variation from seasonal air–sea heat fluxes, advective heat transport from upwelling and tidal circulation, and daily heating and cooling of shallow waters. Despite this, seagrass–temperature studies typically only examine a single isolated temperature process, often seasonal heating/cooling or marine heatwaves. Furthermore, elucidation of relationships between different short-term temperature processes and seagrass metrics could provide insights into biologically relevant temperature metrics for seagrasses. Here, we examine the effects of multiple short-term temperature processes on Zostera marina beds in Atlantic Canada, by describing the seasonal phenology of bed characteristics, plant morphology, and physiology across different temperature regimes and by identifying relationships between different temperature and seagrass metrics. We also include water depth as a proxy for light availability. Four distinct short-term temperature processes (median temperature, growing degree day [heat accumulation], daily temperature range, and time in the optimal temperature range [5–23°C]) were used to categorize temperature regimes across our study sites as warm and highly variable, or cool and less variable. We found that both temperature regime and light availability were important for leaf area index (LAI) and shoot density, which both decreased with increasing depth yet had the strongest seasonal differences and maximum rates of increase (shoot density) or lowest values and maximum rates of decrease (LAI) in warm and highly variable temperature regimes. These temperature regimes were also associated with seasonal patterns in number of leaves, reduced leaf lengths and number of leaves per shoot, and thinner rhizomes with less carbohydrate storage. Generalized additive models (GAMs) for each seagrass metric using the four temperature processes indicated that while median temperature displayed expected relationships with almost every seagrass metric, inclusion of other temperature processes revealed additional insights not evident from median temperature alone. Our study shows that different sources of short-term temperature variation influence seagrass properties and that accounting for these will improve our understanding of the temperature–seagrass interaction.

INTRODUCTION

Seagrass growth, abundance, and distribution are regulated primarily by light availability, water temperature, and nutrients (Lee et al., 2007). While seagrass phenology has been strongly linked to seasonality in these drivers, the importance of temperature in maintaining seagrass populations has been widely recognized, particularly within the context of climate change (Nguyen et al., 2021). Water temperature time series are readily and easily measured and are comprised of many different sources of variability that occur across different timescales (e.g., decadal warming, seasonal warming and cooling, wind-driven upwelling, and tidal temperature fluctuations). In the marine nearshore, localized phenomena such as solar heating, riverine input, meteorological events, and tidal circulation interact with the complex coastline and seabed to regulate temperature on short timescales (hours to months; Wong & Dowd, 2021), resulting in unique temperature signatures at different sites. These temperature regimes can occur across small spatial scales and can have important effects on seagrasses (Krumhansl et al., 2021). To date, most seagrass–temperature studies examine only single short-term temperature processes in isolation, typically focusing on seasonal heating and cooling (represented by mean or median temperatures) or, more recently, marine heatwaves and threshold exceedances (e.g., Marbà & Duarte, 2010; Reusch et al., 2005; Shields et al., 2019; Strydom et al., 2020). Inclusion of multiple short-term temperature processes (e.g., seasonal changes, heat accumulation, and daily range) would not only account for their collective impacts on seagrasses, but also allow specific relationships between each temperature process and various seagrass properties to be identified. This would provide new insights into biologically relevant temperature metrics for seagrasses and help identify potentially important mechanistic pathways worthy of future investigation.

Temperature can influence seagrasses directly. Direct biochemical pathways include photosynthesis and respiration, where photosynthetic rates increase with temperature to certain optimal thresholds (Lee et al., 2007), beyond which respiration will outpace photosynthesis, causing reduced photosynthesis to respiration (P:R) ratios (Bulthuis, 1987; Lee et al., 2007). Physiological responses to thermal stress include enhanced production of carbohydrates, heat shock proteins, and reactive oxygen species, all of which protect cellular components from thermal damage (Gu et al., 2012; Marín-Guirao et al., 2018). Morphological responses include changes to leaf length, number of leaves, biomass allocation, and growth rate, which help optimize carbon balance, while at the population level, seagrasses can either acclimate to local temperature conditions or shift their range distribution (Nguyen et al., 2021). Additionally, temperature can also act in conjunction with other physical variables to affect seagrasses, with light being particularly important. Here temperature strongly modulates the relationship between photosynthetic rate and underwater irradiance, with increased temperature causing increased light saturated photosynthetic rate, compensation irradiance, and saturation irradiance (Bulthuis, 1987; Lee et al., 2007; Marsh et al., 1986). Also, because temperature disproportionately increases respiration relative to photosynthesis, seagrasses exposed to high water temperatures have higher light requirements to maintain carbon balance (Moore et al., 2012).

In addition to direct effects and acting with other variables to affect seagrasses, temperature can also indirectly provide information on other important physical drivers at a site (i.e., as a proxy variable). For example, our previous study showed that deep, exposed, and well flushed sites had cool and less variable temperature regimes compared with shallow, protected, and poorly flushed sites that had warm and highly variable temperature regimes (Krumhansl et al., 2021). These temperature regimes were associated with specific seagrass properties, where seagrass production and resilience were highest at the deep exposed sites (Krumhansl et al., 2021). Here exposure and flushing not only influenced the overall temperature signal but also likely impacted seagrasses through nutrient delivery and sediment properties. Thus, using temperature as a proxy for other physical drivers important for seagrasses allows site-specific conclusions to be drawn that are collectively shaped by the various physical processes present. This is valuable as it enables general in situ conclusions based on temperature regime without having to measure and examine each physical process separately, which are often difficult to disentangle in field settings.

Both long- and short-term temperature processes can influence the dynamics of seagrass ecosystems. Long-term temperature processes that span multiple years or decades and typically have large spatial scales are associated with climate oscillations (e.g., El Niño-Southern Oscillation [ENSO], North Atlantic Oscillation [NAO], and Pacific Decadal Oscillation [PDO]). These long timescales make it difficult to study effects on seagrass ecosystems, given the necessity for coincident long time series data for both temperature and plants. Despite this, processes such as anthropogenically induced ocean warming and the ENSO have been shown to reduce seagrass cover and survival and also cause range shifts and changes to species compositions (Johnson et al., 2003; Marbà & Duarte, 2010; Richardson et al., 2018; Thom et al., 2003; Virnstein & Hall, 2009). In the marine nearshore, these long-term, large-scale temperature trends provide the backdrop upon which short-term, more localized temperature variability is superimposed. Nearshore temperature variability is often dominated by short-term processes that operate on timescales of hours to months (i.e., seasonal and sub-seasonal processes). These processes are highly variable across small spatial scales, owing to interactions between localized phenomena and coastal bathymetry and morphology. Temperature fluctuations at the hourly and daily timescales result from solar heating of shallow waters, low tides that coincide with peak solar insolation hours, tidal exchange between inshore/offshore waters, restricted water flow due to coastal morphology, and tidal influence on flushing properties and currents (George et al., 2018; Krumhansl et al., 2021; Pedersen et al., 2016; Wong & Dowd, 2021). Temperature variations on timescales of days to months result from seasonal heating and cooling as well as wind-driven coastal upwelling (Bylhouwer et al., 2013; Petrie et al., 1987). Although seagrass productivity and resilience are simultaneously affected by several different short-term temperature processes (Krumhansl et al., 2021), few studies have examined the collective impacts of multiple short-term temperature processes on seagrasses (but see Krumhansl et al., 2021; Shields et al., 2019). In fact, most studies focus on a limited set of temperature processes that include either seasonal temperature trends (using mean or median temperature; e.g., Lee et al., 2007; Shields et al., 2019; Wong et al., 2013), marine heat waves, or localized warming (e.g., Collier & Waycott, 2014; Kim et al., 2020; Moore et al., 2014; Richardson et al., 2018).

A full understanding of temperature effects on seagrasses requires simultaneous consideration of the various sources of short-term temperature variation that operate in the nearshore, in addition to other key physical processes such as light availability. Previously we showed that several short-term temperature processes (seasonal trend, tidal and daily temperature range, yearly heat accumulation, and threshold exceedances) along with light availability were primary drivers of seagrass productivity and resilience to disturbance (Krumhansl et al., 2021). However, this work did not include a time series of plant data to examine how short-term temperature variation interacts with light to affect seagrass phenology. Furthermore, it did not identify specific relationships between key short-term temperature processes and seagrass metrics, precluding insights into biologically meaningful temperature metrics for seagrasses. Here, we examine these aspects for Zostera marina (eelgrass) beds on the Atlantic coast of Nova Scotia, Canada. We use time series data of in situ temperature and plant measurements from seagrass beds inhabiting different temperature regimes to (1) describe the phenology of Z. marina bed characteristics, morphology, and physiology across different temperature conditions, (2) characterize key short-term temperature processes evident in the seagrass beds, and (3) determine the relationships of these short-term temperature processes with metrics of eelgrass cover, morphology, and physiology. We then discuss the biological relevance of these temperature processes and how they might inform understanding of the seagrass–temperature interaction and hypothesis generation for focused field and laboratory studies. Note that we include water depth as a proxy for light availability in our analyses, given its known interactive importance with temperature, particularly in shaping the phenological responses of seagrasses.

MATERIALS AND METHODS Study sites

The effects of short-term temperature processes and light availability on seagrass beds were examined at six sites on the Atlantic coast of Nova Scotia, Canada: Port L'Hebert (PH), Lower Three Fathom Harbour (L3F), Sambro (deep and shallow), and Taylor Head (TH; deep and shallow) (Figure 1). These sites were selected as they are subjected to a wide range of physical conditions that shape eelgrass productivity and resilience to disturbance (Krumhansl et al., 2021; Wong & Dowd, 2021). PH and L3F were relatively shallow (1.7 and 0.9 m deep at mean high tide, respectively) and warm, with muddy/silty sediments, and low wind exposure. L3F also had reduced flushing owing to its location within a lagoon with highly restricted tidal flows. In contrast, Sambro and TH were relatively deep and cool, with sandy sediments, high wind exposure, and regular tidal flushing. The seagrass beds at Sambro and TH spanned a larger depth gradient than at PH and L3F, so both a shallow (1.8 and 3.2 m deep at Sambro and TH, respectively) and deep (6.2 and 5.1 m deep at Sambro and TH, respectively) site were sampled at each location. Across all sites, maximum and mean bottom photosynthetically active radiation (PAR) were strongly negatively correlated with depth at mean high tide (r = −0.90 and −0.89, respectively); shallow sites had higher bottom PAR than deeper sites. For further detailed description of the physical conditions at the sites, see Krumhansl et al. (2020, 2021), and Wong and Dowd (2021).

FIGURE 1. Study sites along the Atlantic coast of Nova Scotia, Canada.

Characterizing seagrass phenology

Field sampling was conducted monthly from mid-May 2018 to mid-July 2019, with the exception of January to March 2019 when weather restrictions precluded sampling. On each sampling date, 10 sampling stations within each seagrass bed were haphazardly distributed at approximately the same depth and sampled. Stations were at least 10 m apart and 2 m away from any seagrass-bare sediment interface. Sampling was conducted using snorkeling at PH and L3F and scuba diving at Sambro and TH. Sampling was conducted at approximately the same tide level (mid to low tide) during all sampling periods.

At each sampling station, water depth was recorded and the number of seagrass shoots (vegetative + reproductive) were determined within a 0.25 × 0.25 m quadrat. Three mature vegetative shoots with 5 cm or more live rhizome were collected and refrigerated for 1–2 days prior to processing. Above- and belowground plant biomass was collected from six of the sampling stations using a hand corer (10.8 cm diameter × 12 cm depth) and refrigerated for 2–4 days until processed.

In the laboratory, leaf length was measured as the distance from the insertion point to the leaf tip. Leaf width was measured midway along each leaf. The number of leaves per shoot was determined by counting all leaves (broken and full) that extended above the insertion point, excluding new leaves fully contained within the sheath. At a subset of sites (PH, TH shallow, and TH deep), rhizome width was measured laterally in the middle of the third rhizome internode. A minimum 50-mm healthy and live piece (i.e., light brown and firm with roots removed) was frozen at 80°C for further analysis of water soluble carbohydrates (WSCs) (see below). Plant material collected in hand corers was placed into a mesh bag (10-mm-diameter holes) and rinsed in salt water to remove sediments. Cleaned material was then separated into live aboveground and belowground biomass components, dried at 60°C for 48 h, and weighed, and the ratio of above to belowground biomass determined. Leaf area index (LAI, leaf area [m²] per bottom area [m²]) was determined for each quadrat by multiplying the mean total leaf area (one-sided) per shoot by the shoot density.

Frozen rhizomes were analyzed for WSCs (i.e., soluble sugars). Starch content was not determined since previous analyses and literature indicate that it is very low in Z. marina (<1% of nonstructural carbohydrates in rhizomes; B. Vercaemer, unpublished data; Amber Sadowy, BC Ministry of Environment and Climate Change Strategy). Each rhizome sample was first freeze-dried and ground in a bead mixer mill, and then WSC were extracted by adding 80% ethanol, heating to 90°C for 10 min, and centrifuging at 13,000g for 1 min. WSC in the supernatant was then quantified using the phenol sulfuric acid method in a microplate assay read at 490 nm (Masuko et al., 2005; Wong et al., 2020).

Temperature dynamics

Time series records of water temperature spanning the duration of the study were obtained at each site using temperature data loggers (TidbiT v2, Onset Computer Corporation) deployed on the sea bottom that recorded every 15 min. Temperature metrics that represent distinct short-term processes were calculated from these time series records and included (1) median temperature, (2) growing degree day (GDD), (3) daily temperature range, and (4) the proportion of time spent in the optimal temperature range (OTR) for photosynthesis. All metrics with the exception of GDD were calculated every 3-week period prior to sampling; this accounted for the preceding temperature history that could plausibly mechanistically influence seagrass responses (Nguyen et al., 2020). It also acts to exclude any disproportionate influence of single anomalous events. Temperature loggers at PH and L3F were sometimes exposed to the air at low tide, as is the seagrass bed. This resulted in both higher (>30°C) and lower (below freezing) temperatures than would be expected if the loggers remained submerged. We retained these temperatures because both warm and cold air temperatures have been shown to impact seagrasses (Park et al., 2016).

The temperature metrics were computed as follows. Median temperature was determined as the median of the recorded temperatures for each analysis period (i.e., every 3-week period prior to sampling). Daily temperature range was calculated by first determining the difference between the 90th and 10th percentiles in temperature for each day and then taking the median of these values for the analysis period. Time in the OTR was calculated as the proportion of time spent between 5 and 23°C for each analysis period. This range was chosen because the photosynthesis to respiration ratio (P:R) for Z. marina peaks at 5°C (Marsh et al., 1986) while respiration begins to outpace photosynthesis at 23°C (Lee et al., 2007). GDD was calculated for each sampling year using the approach of Neuheimer and Taggart (2007), where[Image Omitted. See PDF]and[Image Omitted. See PDF]

Here, the GDD at each time interval t (day) is calculated by subtracting a base temperature (T_base) from the mean of the daily maximum (T_max) and minimum (T_min) temperatures, and these values are then cumulatively summed across each yearly temperature record for total heat accumulation. We extracted the accumulated GDD from each sampling date to use in the analyses. We used a base temperature of 5°C, as with the OTR calculations, because this is where the P:R ratio peaks. We did not impose a maximum cap to temperature as is sometimes done for land plants, since the P:R ratio for Z. marina remains above 1 for temperatures up to 35°C (i.e., beyond the maximum temperature observed in this study).

Other metrics describing short-term temperature variability were also considered (e.g., SD, 10th and 90th percentiles, kurtosis, and skew) for summarizing the raw temperature data, as well as the temperature time series separated into different frequency bands, such as the tidal and meteorological bands (Wong & Dowd, 2021). However, these metrics were either highly correlated with other temperature metrics or did not have strong correlations with the plant metrics and so were not used for the analyses. The four temperature metrics retained for our analyses represent relatively distinct metrics of both seasonal (months) and sub-seasonal (hours to weeks) short-term temperature processes. Median temperature captures the gradual seasonal heating and cooling of water. Daily temperature range accounts for variability associated with localized processes, including solar heating, episodic wind-driven upwelling, and tidal flushing. Time in the OTR accounts for variability associated with temperature exceedances outside of the optimal thermal threshold, while GDD represents heat accumulation in a system across the year. In addition to their representativity, we also chose these metrics based on their potential as highly biologically relevant metrics that could be directly linked to seagrass properties and site characterizations and thus provide insights into underlying mechanisms.

Statistical modeling Characterizing temperature regimes

The overall temperature regime at each site was characterized from the magnitude and variability associated with each short-term temperature process examined. Further, we conducted a principal components analysis (PCA) to examine how sites were grouped across different temperature regimes, as defined by the four temperature processes collectively. We also included water depth as a proxy for light availability, given its strong negative correlation with bottom PAR. The largest differences in temperature processes among sites were observed between June and September 2018, and data from these months were used for the PCA. Mean values of each temperature process and depth were calculated within this time period and standardized prior to use in the PCA. The PCA was conducted using the package “vegan” in R Studio (2021.09.1).

Seagrass phenology

We used a generalized additive model (GAM) framework to characterize the phenology of various seagrass metrics and to establish relationships with the different short-term temperature processes and water depth. Additive models such as GAMs allow flexible specification of the relationship between the response and covariate(s) by incorporating a sum of smoothed functions for one or more covariates (Wood, 2017). These models are thus appropriate to describe nonlinear relationships often evident for seagrass metrics over seasons and across environmental conditions (Lefcheck et al., 2017; Wong et al., 2021). Response variables included LAI, total shoot density (vegetative + reproductive shoots), aboveground to belowground (AGBG) biomass ratio, length of leaf 3, number of leaves per shoot, rhizome width, and rhizome carbohydrates.

The phenology of each seagrass metric over the ~15-month sampling period was examined using separate GAM models with site as a fixed categorical predictor and time as the covariate. A gamma distribution with the log link function was used for all metrics except number of leaves per shoot, which used a Poisson distribution with log link, and rhizome WSC, which used a Gaussian distribution with an identity link. The resultant smooth GAM curves of the seagrass metrics over time were used to obtain their maximum and minimum values. Maximum values were determined for the first year of sampling (April 2018 to April 2019). Minimum values were identified for June 2018 to June 2019; this time period allowed identification of winter minimums, as well as minimums that occur outside the winter period (i.e., for the AGBG biomass ratio). Maximum growth and decay rates were also calculated from the first derivative of the local slope of the predicted smooth curves for all plant metrics except for the number of leaves per shoot, which showed only seasonal patterns at two sites.

Relationships of temperature processes and water depth with seagrass metrics

The relationship between each seagrass metric (response variable) and the four different temperature processes and water depth (predictor variables) was examined using generalized additive mixed models (GAMMs). GAMMs are a generalization of GAMs that allow for the incorporation of both fixed and random effects (i.e., mixed models). For each case, site was included as a random effect, and covariates were the temperature metrics (i.e., median temperature, GDD, daily temperature range, and time in the OTR) and water depth, used here as a proxy for light availability. We were unable to use time series of bottom PAR directly as PAR measurements were restricted to the summer period. Anomalies from the common seasonal signal across sites were used for both the response variables and temperature covariates. The use of anomalies helps to isolate differences in temperature variability at each site by removing the confounding common seasonal signal, allowing more interpretable relationships between plant metrics and temperature processes. The anomalies were calculated by fitting a smoothing spline to each plant response variable or temperature covariate over time for data across all sites, then subtracting the predicted common signal from the metric at each sampling date. It follows that positive anomalies indicate the metric is higher than the common signal (or higher than average), while negative anomalies mean the metric is lower than the common signal (or lower than average). We use only the terms “higher than average” and “lower than average” in the text to simplify explanations. The GAMMs all used a Gaussian distribution with identity link function. Because rhizome width and WSC were only examined at three sites (i.e., PH, TH shallow, and TH deep), the models for these response metrics excluded site as a random effect and the covariate OTR, to allow better definition of the remaining covariate relationships. Additionally, GAMMs for morphological metrics, calculated from several shoots collected from each quadrat, did not include quadrat as a random effect nested within site, as initial analyses indicated that quadrat did not explain a large amount of variance in the data.

All GAM(M) regression analyses for seagrass phenology and relationships with temperature metrics and water depth were conducted using the “gam” implementation in the package “mgcv” in R Studio (2021.09.1). For all regressions, the smoothing terms relied on a shrinkage version of the cubic regression spline, which reduces the coefficient of non-influential covariates towards zero (Wood, 2017). The random effect site was included in the GAMM models using a smoother function (Wood, 2017). For all models, the significance of smoothers (set as p < 0.05) was evaluated using p-values from the regression statistics. Model fit was assessed using the percent deviance explained. Residual plots were evaluated to assess the assumptions of homogeneity of variance, normality, and temporal autocorrelation; no violations were evident.

RESULTS Temperature conditions across sites

The PCA indicated that sites were strongly separated based on their overall temperature regime and water depth (Figure 2). The first principal component (PC1) explained the majority of the variance across sites (90%) and was defined using all four temperature processes and water depth. The small proportion of variance explained by the remaining PCs suggests the temperature metrics used along with water depth provide a fairly compact and complete representation of site conditions. Sites on the negative side of PC1 (L3F and PH) were defined by higher daily temperature ranges, median temperature, and heat accumulation (GDDs), lower time spent in the OTR, and lower water depth. Sites on the positive side of PC1 (shallow and deep sites at Sambro and TH) had lower daily temperature ranges, median temperature, GDD, more time in the OTR, and higher water depth. Thus, we characterized PH and L3F as sites with warm and highly variable temperature regimes, while Sambro and TH (deep and shallow) were sites with cool and less variable temperature regimes. Water depth indicates that the warm sites had higher bottom light relative to the cooler sites during the summer period.

FIGURE 2. Principal components analysis of the four short-term temperature processes and depth characterized at each site. Loadings of temperature processes and water depth are indicated by the ends of the vectors. OTR, optimal temperature range; PC, principal component; Prop., proportion.

This classification was further supported by evaluation of the individual temperature processes. Median temperatures at PH and L3F during the growing season were consistently higher than other sites by 5–10°C (Figure 3a) and often exceeded 20°C. Heat accumulation was consistently highest at PH and L3F, with heat accumulating earlier in the season, accumulating at higher rates throughout the year (i.e., higher slopes), and reaching higher values in the late fall (Figure 3c). Daily temperature ranges were also higher at PH and L3F than other sites throughout the year, being 3–4°C higher in the summer and 1–2°C in the winter (Figure 3b). Finally, PH and L3F spent only 40%–80% of the time within the OTR (5–23°C) during the growing season (Figure 3d), with remaining time spent above the upper optimal threshold. At Sambro and TH, the shallow sites were consistently warmer than deeper sites, had higher heat accumulation, and sometimes had higher daily temperature ranges. However, these differences between shallow and deep sites were much smaller than for comparisons with L3F and PH. Thus, the combination of these observations with the PCA results supports our classification of PH and L3F as sites with warm and highly variable temperature regimes, and Sambro and TH as sites with cool and less variable regimes.

FIGURE 3. The four short-term temperature processes examined, including (a) median temperature, (b) daily temperature range, (c) growing degree day (GDD), and (d) proportion of time in the optimal temperature range (OTR). All processes were calculated for the 3-week period prior to each sampling event, with the exception of GDD, which was calculated as a cumulative sum across each year. The green line is the smoother fit to the data, from which anomalies were calculated. Tick marks along the x-axes indicate sampling dates. See text for details. Prop., proportion.

Phenology of seagrass Seagrass bed characteristics

Strong seasonal patterns in LAI were evident across all sites, regardless of temperature regime or depth (site smooth terms: p < 0.001, model deviance = 73.2%; Figure 4a). Maximum LAI was similar across all shallow sites regardless of temperature regime (5.9–7.2 m² m⁻²; Appendix S1: Table S1), although it occurred earlier at sites where temperatures were warm and highly variable (i.e., PH and L3F). Maximum LAI was also consistently higher at all shallow sites compared to deep sites. In contrast, lowest LAIs were observed only at shallow sites with warm and highly variable temperature regimes (i.e., 0.40 and 0.53 m² m⁻² at L3F and PH, respectively). The maximum rate of increase was highest at PH, followed by the remaining shallow sites. PH and L3F exhibited the highest maximum rates of decrease in LAI that occurred earliest in the year. Interestingly, observation of the phenology in shoot density and leaf length (Appendix S1: Figures S1 and S2) indicates that phenology of LAI at shallow sites (both warm and cool) was driven by seasonal changes in both shoot density and leaf length. However, at deeper sites (TH deep, Sambro deep), phenology of LAI was driven primarily by seasonality in leaf length, as shoot density did not have seasonal patterns at these sites.

FIGURE 4. Phenology of (a) leaf area index (LAI) and (b) aboveground to belowground (AGBG) biomass ratio at each field site estimated by the generalized additive model (black line) using the monthly plant data (colored circles). Gray shading indicates CIs based on 2 SEs.

The strong seasonal patterns in shoot density observed at shallow sites coupled with their absence at deep cool sites is an important contrast to the phenology of LAI (site smooth terms: p < 0.006, model deviance = 62.1%; Appendix S1: Figure S1). The largest difference between minimum and maximum shoot density and also the highest rate of increase (4.3 shoots day⁻¹) was observed at PH, where temperature regime was both warm and highly variable and water depth was shallow (Appendix S1: Table S1). At Sambro and TH shallow, where temperatures were cooler and less variable and water depths deeper than PH or L3F, maximum rates of increase were much lower than at PH but occurred 1–2 months earlier (1.1 and 0.39 shoots day⁻¹, respectively). Although L3F exhibited strong seasonal phenology, maximum growth rate is not reported as it was not properly captured in the first sampling year. Maximum rates of decrease in shoot density were earlier at sites with warmer and more highly variable temperature regimes (PH and L3F) compared with cooler and less variable sites (Sambro and TH, both deep and shallow).

Similar to LAI, AGBG biomass ratio showed seasonal cycles across all sites, regardless of depth or temperature. However, the strongest seasonal cycles were observed at sites with warm and highly variable temperature regimes (L3F and PH) compared with sites with cooler and less variable regimes (Sambro and TH), regardless of depth (site smooth terms: p < 0.01, model deviance = 51.5%; Figure 4b). AGBG biomass ratio values at L3F and PH peaked earlier in August to September compared with other sites, before declining throughout the fall and winter (Appendix S1: Table S1). While seasonal cycles in AGBG biomass ratio were evident at Sambro and TH, differences between maximum and minimum values were less, and the cycle amplitudes were dampened. Maximum rates of both increase and decrease were highest at PH and L3F relative to Sambro and TH.

Seagrass morphology and physiology

Leaf length at all sites was strongly seasonal regardless of temperature regime or water depth, with lengths increasing from spring to early fall before declining during late fall and winter (site smooth terms: p < 0.001, model deviance = 76.3%; Appendix S1: Figure S2). Maximum leaf length occurred earliest at L3F in August, relative to all other sites where leaf lengths were longest in October to November (Appendix S1: Table S1). Maximum leaf length tended to be shorter at the shallow sites with warm and highly variable temperature regimes, although it was similar also at TH shallow.

Seasonal trends in the number of leaves per shoot were only observed at shallow sites with warm and highly variable temperature regimes (L3F smoother term: p < 0.001, PH smoother term: p = 0.002, model deviance = 28.1%; Figure 5). A strong decline in the number of leaves across the entire sampling period was observed at Sambro deep (smoother term: p = 0.005), while it did not change across the sampling period at the remaining cool sites (smoother terms: p > 0.05). Highest and lowest number of leaves per shoot were both observed at L3F and PH.

FIGURE 5. Phenology of number of leaves per shoot at each field site estimated by the generalized additive model (black line) using the monthly plant data (colored circles). Gray shading indicates CIs based on 2 SEs.

Rhizome width was measured at a subset of sites and showed seasonal patterns at both PH and TH shallows, while it remained constant at TH deep (PH and TH smoother terms: p < 0.001, model deviance = 70.3%; Figure 6a). The temperature regime appeared to affect the patterns observed at the shallow sites. At TH shallow, where the temperature regime was cool and less variable, rhizome width increased during the spring and into the fall before declining over the winter (Appendix S1: Table S1). At PH, where temperature regime was warm and highly variable, rhizome width declined from the spring into fall and the following winter, increasing only in the following spring. Rhizomes also remained thinnest at PH relative to TH throughout the year.

FIGURE 6. Phenology of (a) rhizome width and (b) rhizome water soluble carbohydrates (WSC) at three selected field sites estimated by the generalized additive model (black line) using the monthly plant data (colored circles). Gray shading indicates CIs based on 2 SEs. DM, dry mass.

Rhizome WSCs showed strong seasonal patterns at all three sampled sites regardless of temperature regime (site smoother terms: p < 0.005, model deviance = 25.4%, Figure 6b). However, rhizome WSC remained constant from April to October 2018 at PH, where temperature regime was the warmest with highest variability, whereas it increased at the cooler sites (TH deep and shallow) during the same time period. Peak rhizome WSC concentrations and maximum rates of increase were highest at TH relative to PH (Appendix S1: Table S1). At all sites, rhizome WSC declined across the winter and then increased the following spring and summer.

Relationships of seagrass metrics with short-term temperature processes and water depth Seagrass bed characteristics

GAMM regressions highlighted the importance of both water depth and short-term temperature processes on patterns in LAI (model deviance = 33.8%; Figure 7). LAI decreased with water depth, being highest at the two most shallow warm sites (L3F and PH) and lowest at the deepest cooler sites (TH and Sambro deep; smoother term: p < 0.0001). LAI increased moderately as median temperature increased across the depth range of the sites (smoother term: p = 0.026), likely illustrating the modulating effect of temperature on light availability. The direct effects of temperature on LAI are most evident from relationships with GDD and daily temperature range. Here LAI exhibits an accelerated decline when certain thresholds in these temperature processes are exceeded (i.e., anomalies of 500 GDD and 1 daily temperature range; smoother terms: p < 0.0001) at the shallow warm sites. This occurs despite LAI remaining relatively constant across the depth range of these sites. Further, LAI also tended to decrease as time in the OTR increased (smoother term: p < 0.0005), with highest values driven by plants at the shallow warm sites.

FIGURE 7. Generalized additive mixed regression model (black lines) illustrating the relationship of leaf area index (LAI) anomalies (colored circles) with anomalies of each short-term temperature process and water depth. Gray shading indicates CIs based on 2 SEs.

Similar to LAI, shoot density also declined across water depth (smoother term: p < 0.0001; model deviance = 45.1%; Appendix S1: Figure S3), suggesting the structuring influence of light availability. However, all temperature processes, with the exception of daily temperature range, were also important for shoot density (smoother terms: p < 0.0001). Shoot density had a convex relationship with median temperature, where it first increased at lower temperature anomalies and then decreased at higher anomalies. Decreased shoot density was driven by plant responses at the warm shallow sites (PH and L3F), which had relatively constant density across the depths sampled. Shoot density was relatively stable at the low heat accumulation values observed at Sambro and TH but then increased as heat accumulated at PH and L3F before declining when GDD anomaly exceeded 500. Shoot density also increased well above average when the proportion of time in the OTR decreased below 40% (smoother term: p < 0.0001). Daily temperature range did not have a significant relationship with shoot density (smoother term: p = 0.382).

In contrast to both LAI and shoot density, the AGBG biomass ratio was influenced only by the short-term temperature processes and did not change across depth (model deviance: 33.8%; Figure 8). The AGBG ratio increased as median temperature increased, with highest ratios evident at PH and L3F (smoother term: p = 0.0008). The AGBG ratio also increased when daily temperature range was positive, but only to a certain threshold (i.e., 1°C above average), above which it declined rapidly, primarily at PH and L3F. The AGBG ratio was not significantly influenced by GDD or OTR (smoother term: p > 0.05).

FIGURE 8. Generalized additive mixed regression model (black lines) illustrating the relationship of anomalies of the aboveground to belowground (AGBG) biomass ratio (colored circles) with anomalies of each short-term temperature process and water depth. Gray shading indicates CIs based on 2 SEs.

Seagrass morphology and physiology

As expected, leaf length was strongly influenced by water depth, increasing as depth increased (smoother term: p < 0.0001; model deviance = 46.2%; Appendix S1: Figure S4). However, water temperature also played a role, evident from plant responses at PH and L3F to all temperature processes except to daily temperature range (smoother terms: p < 0.0001). At average to slightly below-average median temperatures, leaf lengths were relatively stable. However, they increased dramatically when median temperature was 0–4°C above average, after which it decreased rapidly, a pattern driven by PH and L3F. Similarly, when heat accumulation exceeded 500, leaf lengths at PH and L3F decreased rapidly to below average values. Leaf lengths were relatively stable across the range of OTR at all sites, except for at the highest OTR values at Sambro and TH where leaves were longer than average (smoother term: p < 0.0001). No relationship was evident between daily temperature range and leaf length (smoother term: p = 0.128).

The GAMM for the number of leaves per shoot did not explain a large portion of the variance in the data (model deviance = 19.2%); however, it still provides insight into relationships with temperature processes (Figure 9). No pattern in number of leaves was evident across the depth range sampled (smoother term: p = 0.194). However, the number of leaves tended to decrease at higher median temperatures, with more heat accumulation, and when less time was spent in the OTR (smoother terms: p < 0.0001). Daily temperature range did not influence the number of leaves (smoother term: p = 0.214).

FIGURE 9. Generalized additive mixed regression model (black lines) illustrating the relationship of anomalies of the number of leaves per shoot (colored circles) with anomalies of each short-term temperature process and water depth. Gray shading indicates CIs based on 2 SEs.

The GAMM for rhizome width indicated that only heat accumulation and daily temperature range were important (smoother terms: p < 0.0001, model deviance = 57.9%), although the data suggest rhizome width also increased as depth increased (smoother term: p = 0.298; Figure 10). Rhizome width had a convex relationship with GDD, with rhizomes being thinner than average at PH when heat accumulation was very high (GDD anomaly = 750). Rhizome widths also tended to decrease as daily temperature range increased at PH, with rhizomes being thinner than average when daily temperature range was high. Rhizome widths did not vary with median temperature or depth (smoother terms: p = 0.516 and 0.298, respectively).

FIGURE 10. Generalized additive mixed regression model (black lines) illustrating the relationship of rhizome width anomalies (colored circles) with anomalies of each short-term temperature process and water depth. Gray shading indicates CIs based on 2 SEs.

The GAMM for rhizome WSCs did not explain a large portion of the variance in the data (model deviance = 26.2%). However, the smoothers for median temperature and daily temperature range provide some insight into relationships with rhizome WSC (smoother terms: p < 0.0001 and p = 0.0008, respectively; Figure 11). The peak in rhizome WSC in the median temperature smoother is likely an artifact of missing data (i.e., interpolation of a smooth function over a gap). However, WSC did initially increase at PH with positive median temperature anomalies to a threshold (5°C), after which it declined. Rhizome WSC tended to increase as the daily temperature range increased from 0 to 1°C above average but then subsequently decreased above this value. GDD and depth did not have significant relationships with rhizome WSC (smoother term: p = 0.06 and 0.397, respectively).

FIGURE 11. Generalized additive mixed regression model (black lines) illustrating the relationship of rhizome water soluble carbohydrate anomalies (colored circles) with anomalies of each short-term temperature process and water depth. Gray shading indicates CIs based on 2 SEs. DM, dry mass.

DISCUSSION

We showed that the seagrass beds studied inhabited temperature regimes that ranged from relatively cool temperatures with low variability to relatively warm temperatures with high variability. These temperature conditions explained patterns in the phenology and overall properties of the seagrass beds. We focused on both temperature and light (using depth as a proxy) because they are considered primary drivers for seagrasses (Krumhansl et al., 2021; Lee et al., 2007). Both variables were particularly important for bed characteristics and leaf lengths, which reflected the seasonal patterns of light and temperature. While it can be difficult to disentangle the individual role of temperature and light (Lee et al., 2007), our observations across the different temperature regimes coupled with the specific temperature and depth relationships with different seagrass metrics highlighted the potential role of temperature in shaping seagrass ecosystems. Furthermore, we showed that while median temperature was important for almost every seagrass metric explored, heat accumulation, daily temperature range, and time in the OTR provided additional insights. This suggests their potential value as biologically relevant temperature metrics for seagrasses. We further explore these findings below.

The phenology and overall properties of seagrass bed characteristics, particularly LAI and shoot density, were related to both temperature regime and light availability. The role of light availability was evident when observing phenological patterns across sites coupled with the relationships of LAI and shoot density with depth identified in the GAMM analyses. Although strong seasonal patterns in LAI were observed across all temperature regimes and depths, maximum LAI was consistently highest at shallow sites, regardless of temperature conditions. Examination of shoot density and leaf length phenology, the main components of LAI, indicated that both were important for LAI at shallow sites, while leaf length was the primary driver at the deeper sites. Although our observed negative relationships of shoot density and leaf length with depth have been commonly observed, the nonlinear decline in LAI with depth has been less frequently characterized (but see Enríquez & Pantoja-Reyes, 2005; Enríquez et al., 2019; Ralph et al., 2007). LAI is a valuable measure to interpret potential seagrass bed productivity, given it represents the amount of leaf tissue available for light capture (Enríquez et al., 2019). Lower LAI occurs at higher depths to allow an open canopy that facilitates light penetration in order to maintain photosynthetic activity under lower light conditions, while higher LAI in shallow depths is thought to protect plants from photo-inhibition and UV damage (Enríquez et al., 2019). These relationships were reflected in our study and illustrate the important role of light availability in structuring seagrass beds.

We further found that seagrass bed characteristics were also shaped by temperature. Minimum LAI and its maximum rate of decrease, maximum differences in shoot density and its maximum rate of increase, and strongest seasonal cycles in the AGBG biomass ratio and its maximum rates of increase and decrease were all associated with warm and highly variable temperature regimes. Furthermore, LAI, shoot density, and AGBG biomass ratio at warm sites all remained relatively stable over depth but had strong responses to the various temperature processes. The simultaneous evaluation of different temperature processes allowed insights not evident from the most basic and commonly used temperature metric (i.e., median temperature) alone. For example, the convex relationship between median temperature and shoot density illustrates the temperature dependence of photosynthesis, where Z. marina growth, shoot production, and survival are reduced when temperature thresholds are exceeded (Marsh et al., 1986; Nejrup & Pedersen, 2008). These reductions are further reflected at high heat accumulation. However, although we observed reductions in shoot density at high median temperatures and GDD, rapid increases in shoot density were actually evident with intermediate heat accumulation and when temperatures exceeded the optimal threshold 30% or more of the time. This enhanced growth response to anomalously warm temperatures likely represents a long-term shift to greater clonal reproduction rather than maintenance of parent shoots (DuBois et al., 2020; Krumhansl et al., 2021). Similarly, understanding of the temperature–LAI and AGBG biomass ratio relationships was also enhanced by examining multiple sources of temperature variation. LAI increased with median temperature but declined sharply at high heat accumulation and daily temperature ranges, indicating a similar response to high temperature as for shoot density. The AGBG biomass ratio increased with median temperature, suggesting reduced production of belowground tissues to lessen respiratory burden but then further declined at high daily temperature ranges. This indicates that both above- and belowground tissue production can be impacted by large temperature fluctuations (Lee et al., 2005, 2007), a feature not evident from examining median temperature alone.

Accounting for multiple short-term temperature processes also enhanced understanding of shoot leaf structure and rhizome morphology and physiology. Seasonal patterns in the number of leaves per shoot were only evident in warm and highly variable temperature regimes, where leaf loss occurred with high median temperatures and high heat accumulation, but leaf gain occurred when 40% of the time or less was spent outside the OTR. Leaf lengths, which increased with depth, showed gains in length at intermediate and positive median temperature anomalies, while strong reductions in length were evident at high median temperatures and heat accumulation. Reductions in leaf number and length at high temperatures is a common mechanism employed by thermally stressed seagrass to reduce carbon requirements, while increases are often observed at intermediate temperatures (Lee et al., 2007; York et al., 2013). Rhizomes were thinner overall and decreased in width across the summer at warm sites, with significant declines at high daily temperature ranges (1.5–2°C anomaly) and high heat accumulation (>250 GDD anomaly). This suggests reallocation of tissue production to reduce belowground respiratory burden (Lee et al., 2005), although thinner rhizomes typically transport less resources (Duarte, 1991). WSC concentrations at the warm site were stable over the summer, suggesting that excess carbohydrates produced were not stored but used to compensate for high carbon requirements in warm conditions (Wong et al., 2020). This differed from the cooler site where excess WSC were produced and stored over the summer, as is common for temperate non thermal stressed Z. marina (Burke et al., 1996; Vichkovitten et al., 2007; Wong et al., 2020). Relationships of WSC with median temperature indicated reduced production and/or mobilization at intermediate temperatures, but enhanced production at higher median temperatures—a response likely invoked to protect cellular components (Gu et al., 2012; Krumhansl et al., 2021; Moreno-Marín et al., 2018). However, this enhanced production is restricted at anomalously high daily temperature ranges. The response of rhizomes to temperature likely contributed to phenological patterns in the AGBG biomass ratio, which showed the strongest seasonal cycles at sites with warm and highly variable temperature regimes. Given that the belowground component is essential for seagrass resilience to disturbance (Collier et al., 2020), our results suggest that bed resilience will depend in part on the temperature dynamics seagrass beds inhabit.

While we focused on the primary role of temperature and light in shaping seagrass bed properties, other factors such as nutrients and sediment conditions may have also been important. In our study, consistent summer nitrogen content in leaves (Krumhansl et al., 2021) and minimal anthropogenic nitrogen loading (Murphy et al., 2022; Nagel et al., 2018) across our sites suggest that differences in nutrient availability did not influence our results. However, sediments at PH were significantly higher in organic matter content than at remaining sites (20.3% sand and 19.7% organic matter at PH vs. 72.1%–87.9% sand and 1.7%–2.1% organic matter at other sites; Krumhansl et al., 2021), suggesting potential chronic sulphide exposure at PH. Sulphide toxicity in seagrasses can reduce photosynthetic rates, shoot density, root growth, and biomass or even cause extensive die-offs (Borum et al., 2005; Goodman et al., 1995; Pedersen & Kristensen, 2015). Although seagrasses can reduce sulfide intrusion using oxidized micro-zones around the rhizosphere and reduce toxicity in tissues through internal oxidization (Hasler-Sheetal & Holmer, 2015), these mechanisms are strongly influenced by high water temperatures, which can reduce water column dissolved oxygen, enhance sediment sulfide production, or cause physiological stress in the plants (Koch & Erskine, 2001; Pedersen et al., 2004). In our study, thin rhizomes and low belowground biomass observed at PH may have resulted from the effects of both high temperatures and sediment organic content. Reductions in leaf number and lengths at anomalously high temperatures may have further reduced plant ability to oxidize sulfides in both sediments and tissues. Here, temperature likely affected not only sediment processes but also the plant responses.

The temperature regimes observed in our study resulted from multiple temperature processes that were shaped by localized differences in depth, exposure, and water movement. For example, the warm and highly variable temperature conditions at PH and L3F occurred because of shallow water depths which were readily solar-heated during the day and cooled at night, low wind exposure that restricted horizontal mixing and allowed for spatial gradients in temperature and reduced tidal advection that restricted inward transport of cooler offshore waters in the summer. Here these physical processes not only affected the temperature signal but would have also had their own effects on seagrass performance, such as reduced nutrient delivery in low flow environments. We contend that temperature regime is a useful proxy to represent overall site conditions and that while many physical processes shape the temperature signal, it remains associated with specific seagrass properties. Thus, this approach is valuable in that it allows general and biologically relevant conclusions on seagrass responses to different temperature regimes without having to disentangle multiple physical processes, which is not easily achieved for in situ studies.

We previously showed that several temperature processes and bottom light availability were primary drivers of Z. marina productivity and resilience at peak summer biomass, with other physical variables (i.e., sediment type, bottom slope, exposure, and water currents) being of secondary importance (Krumhansl et al., 2021). Here we extended that work by focusing only on temperature and light, using time series plant and temperature data to examine temperature associations with seagrass phenology, and by examining specific relationships between seagrass metrics and individual temperature processes. This approach has seldom been used in previous seagrass–temperature studies but proved powerful in highlighting potentially biologically relevant temperature processes for seagrasses beyond mean/median temperature. For example, LAI strongly decreased with high heat accumulation and daily temperature range, suggesting that these processes deserve further investigation to reveal underlying mechanisms and to identify threshold values that could be useful within a management context. Our choice of temperature processes was based on representing different aspects of its temporal variability associated with seasonal trends, yearly heat accumulation, daily ranges, and thermal exceedances. The interpretability of our analyses was strengthened by using anomalies of both temperature and plant metrics. These removed the common or “average” signal across sites, isolating differences in temperature processes and allowing for more clear and direct interpretations not confounded by seasonal signals common across all sites. Our analyses were further strengthened by characterizing the temperature processes based on their preceding three-week temperature history, which allows temperature metrics to potentially act as casual drivers. Temperature history is important for Z. marina, as it can develop stress memory in response to stressful conditions (Nguyen et al., 2020), allowing epigenetic modifications that improve plant responses to reoccurring disturbances.

In summary, our study showed that Z. marina bed characteristics, morphology, and physiology are strongly affected not only by light availability but also by multiple sources of seasonal and sub-seasonal temperature variability. These temperature processes act across both short timescales and small spatial scales in the marine nearshore, creating heterogenous temperature landscapes that influence seagrass properties, with likely subsequent effects on overall growth and productivity. Furthermore, these short-term temperature processes will shape the local manifestation of larger scale processes such as ocean warming and marine heat waves that are considered key climate change threats to seagrasses (Duarte et al., 2018). The results of our study are likely most relevant for Z. marina in the central part of their distribution, as seagrasses from widely different latitudes can have different molecular responses to environmental stress (Jueterbock et al., 2021; Winters et al., 2011). Nevertheless, our study shows that simultaneous consideration of multiple short-term temperature processes along with light availability can enhance our understanding of seagrass–temperature interactions and provide insights into biologically relevant temperature metrics for seagrasses that should be further explored in the future. Such insights can inform optimal management and conservation strategies that account for the heterogenous nature of nearshore conditions and its response to large-scale environmental change.

ACKNOWLEDGMENTS

We thank B. Vercaemer, J. O'Brien, G. Griffiths, S. Roach, B. Roethlisberger, and M. Scarrow for field and laboratory support. We are grateful for comments from reviewers that improved the manuscript. Funding was provided by Fisheries and Oceans Canada (Melisa C. Wong) and a Natural Sciences and Engineering Research Council of Canada Discovery Grant (Michael Dowd).

CONFLICT OF INTEREST STATEMENT

The authors declare no conflicts of interest.

DATA AVAILABILITY STATEMENT

Data (Wong, 2023) are available on the Government of Canada Open Data Portal: https://open.canada.ca/data/en/dataset/bf0b87d3-f172-4a82-8bc7-4c1a2e0f3ab4.