Adaptations by Zostera marina Dominated Seagrass Meadows in Response to Water Quality and Climate Forcing

: Global assessments of seagrass declines have documented accelerating rates of loss due to anthropogenic sediment and nutrient loadings, resulting in poor water quality. More recently, global temperature increases have emerged as additional major stressors. Seagrass changes in the Chesapeake Bay, USA provide important examples of not only the effects of human disturbance and climate forcing on seagrass loss, but also meadow recovery and resiliency. In the York River sub-tributary of the Chesapeake Bay, the meadows have been monitored intensively using annual aerial imagery, monthly transect surveys, and continuous water quality measurements. Here, Zostera marina has been demonstrating a shift in its historical growth patterns, with its biomass peaking earlier in the growing season and summer declines beginning earlier. We found an increasing trend in the length of the most stressful high temperature summer period, increasing by 22 days since 1950. Over the past 20 years, Z. marina’s abundance has exhibited periods of decline followed by recovery, with recovery years associated with greater spring water clarity and less time spent at water temperatures > 28 ◦ C. Although human disturbance and climatic factors have been altering these seagrass meadows, resilience has been evident by an increase in reproductive output and regrowth from Z. marina seedlings following declines, as well as expansions of Ruppia maritima into areas previously dominated by Z. marina.


Introduction
Recent global climate changes have affected a wide variety of plant communities with diverse geographical distributions. Climate effects on land-based communities have been well studied as to their changing phenologies, ranges, species distributions, and community interactions [1][2][3][4], as well as the additional effects of episodic disturbance factors [5]. Marine seagrass communities have, however, received much less attention until more recently.
Seagrass meadows have been identified as being among the most threatened ecosystems across the globe, as they have been experiencing not only declines in areal extent, but also accelerated rates of decline [6,7]. Both direct and indirect stressors have been implicated in these losses, including increased anthropogenic sediment and nutrient loadings resulting in poor water quality conditions [6,8,9]. Most recently, there has been an increased focus on the effects of current and future climate change factors on seagrass productivity and distribution [10][11][12].
In the Chesapeake Bay, USA, large-scale seagrass declines, first documented in the 1960s and continuing into the 1970s, were related to excessive nutrient and sediment loadings into the Bay, augmented by a devastating tropical storm in 1972 [13]. These declines generated great interest in restoring these critical habitats, and specific targets for water quality and water clarity improvements were put in place to promote seagrass recovery [14][15][16]. The seagrass Zostera marina (eelgrass) has been the dominant species in the polyhaline regions of the Bay, and in the years subsequent to 1972, it disappeared from over 50% of its former Bay-wide distribution [12,13]. In the decades since the 1970s there was some initial recovery [17], however over the past 20 years meadow distribution and abundance have generally been declining, characterized by periods of sharp declines followed by periods of some recovery [16]. Since Z. marina meadows in the Chesapeake Bay are growing near their southern limits of distribution along the eastern coast of North America [18], declines and recoveries in this region can serve as a future model for other areas throughout the northern hemisphere where the interactions of climate and anthropogenic factors, although also very significant, may be less extreme [19,20].
Mean summer water temperatures have been increasing in the Chesapeake Bay, along with greater frequencies of extreme temperature events exceeding 28 • C [12], which has been identified as a critical threshold for Z. marina [21,22]. In 2005 and 2010, extreme heat events triggered large-scale Z. marina die-offs throughout the Bay [20,23]. Z. marina can also be affected by the interactions of temperature and light, with more light required at warmer temperatures to maintain a positive carbon balance [24,25]. Additionally, high temperatures and increased turbidity can reduce oxygen production in the plant, resulting in anoxic conditions which can lead to sulfide toxicity [26,27]. Therefore, the combined effects of climate-related temperature increases and decreased water clarity associated with coastal development will be critical factors in the persistence of this species, not only in the Chesapeake Bay, but throughout its global distribution.
One seagrass-dominated region of the Chesapeake Bay, the York River, has been extensively studied for both water quality conditions and seagrass growth and abundance changes over the past 30+ years, allowing us to combine long term water quality measurements with current (2004-2017) and historical (1980s and 1990s) seagrass abundance surveys to assess Z. marina adaptation to changing conditions. In this study, we used both monthly seagrass transect surveys and annual aerial mapping surveys, combined with both continuous and bimonthly water quality measurements in the York River, Chesapeake Bay to address the following questions: (1) In the face of a changing climate, has there been a shift in the temporal and spatial growth patterns of Z. marina? (2) Can we identify critical climate and water quality conditions associated with periods of Z. marina meadow decline compared with periods of recovery? (3) Although periodic die-off events have been occurring, are these meadows exhibiting signs of or mechanisms for resiliency?

Seagrass Transect Surveys
Seagrass transect surveys were conducted at Goodwin Islands (37 • 13 1 N, 76 • 23 19 W) in the polyhaline portion of the York River Estuary, a tributary of the Chesapeake Bay ( Figure 1). The Goodwin Islands are a 315 ha archipelago of salt-marsh islands that experience semi-diurnal tides with an average range of 0.7 m. The subtidal flats surrounding the islands have supported seagrass beds for at least the last 80 years [13,17,28]. Three permanent transects (GI-1, GI-2, GI-3) were established within this reserve area of the Chesapeake Bay National Estuarine Research Reserve in Virginia (CBNERR) in 2004 as part of their System-Wide Monitoring Program (SWMP). Over this period, divers obtained vegetative cover estimates, maximum shoot lengths, and vegetative and reproductive shoot density counts at monthly intervals along the transects from April-October. The longest transect (GI-3, 700 m) was sampled every 20 m, and the other transects (GI-1, 130 m; GI-2, 300 m) every 10 m. The depth along the transects ranges from 0 to −120 cm below mean lower low water (MLLW). Methods for this survey have been previously published [20]. Biomass was calculated based on shoot length to shoot weight relationships previously developed for this area [29] (shoot weight = shoot length × 0.0018), and then scaled to meter squared estimates based on field density measurements (biomass g DW m −2 = shoot weight × density m −2 ). Z. marina biomass estimates are reported for the entire time period, while percent cover estimates comparing Z. marina and Ruppia maritima (widgeongrass) are reported beginning in 2010, since this is the year that consistent R. maritima data collection began. Similar data from the same location were obtained from historical surveys from the 1980s and 1990s for comparisons [17,29].
Diversity 2018, 10, x FOR PEER REVIEW 3 of 15 maritima (widgeongrass) are reported beginning in 2010, since this is the year that consistent R. maritima data collection began. Similar data from the same location were obtained from historical surveys from the 1980s and 1990s for comparisons [17,29].

Aerial Surveys
Aerial imagery obtained through the Virginia Institute of Marine Science (VIMS) Submersed Aquatic Vegetation Monitoring Program [30] was utilized to quantify seagrass bed area for the Mobjack Bay segment, which includes the Goodwin Islands ( Figure 1). Seagrass bed mapping has been conducted through this program since 1984, using panchromatic photography at a scale of 1:24,000 acquired with a standard mapping camera until 2014. Since 2014, a digital mapping camera with a ground sample distance of 24 cm has been used. The imagery is obtained at peak growing season, typically in June for the polyhaline segments. Guidelines for its acquisition include factors related to tidal stage, plant growth, sun angle, atmospheric transparency, water turbidity, and wind. Complete details of the survey can be found on the website [30].
For the water quality analyses, we limited the seagrass coverage data to 1998-2017, because these are the years water quality data were available. Two distinct periods of decline and two distinct periods of increase were determined from this aerial coverage data. Periods where seagrass was declining included the years 2001-2006 and 2009-2012, while periods of increase included the years 2006-2009 and 2012-2015. These data represent beds during peak Z. marina biomass since they are flown in June, and not peak R. maritima biomass, which is typically low in June and at a maximum in August [17,20].

Water Quality
A long-term temperature record was assembled using samples collected from the Virginia Institute of Marine Science Pier located at 37.2468, −76.5000  and Yorktown USCG

Aerial Surveys
Aerial imagery obtained through the Virginia Institute of Marine Science (VIMS) Submersed Aquatic Vegetation Monitoring Program [30] was utilized to quantify seagrass bed area for the Mobjack Bay segment, which includes the Goodwin Islands ( Figure 1). Seagrass bed mapping has been conducted through this program since 1984, using panchromatic photography at a scale of 1:24,000 acquired with a standard mapping camera until 2014. Since 2014, a digital mapping camera with a ground sample distance of 24 cm has been used. The imagery is obtained at peak growing season, typically in June for the polyhaline segments. Guidelines for its acquisition include factors related to tidal stage, plant growth, sun angle, atmospheric transparency, water turbidity, and wind. Complete details of the survey can be found on the website [30].
For the water quality analyses, we limited the seagrass coverage data to 1998-2017, because these are the years water quality data were available. Two distinct periods of decline and two distinct periods of increase were determined from this aerial coverage data. Periods where seagrass was declining These data represent beds during peak Z. marina biomass since they are flown in June, and not peak R. maritima biomass, which is typically low in June and at a maximum in August [17,20].

Water Quality
A long-term temperature record was assembled using samples collected from the Virginia Institute of Marine Science Pier located at 37.2468, −76.5000  and Yorktown USCG Training Center  (2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017). Sample frequency at the two stations ranged from 2 samples a day to one sample every 6 minutes. The summer growing season for Z. marina was defined based on work conducted in the York River using in situ growth experiments that related net plant growth to water temperature [14] (pp. 95-96). In that study, periods of high growth of Z. marina were observed during the spring and fall, and the low growth summer season was triggered when water temperatures reached 23 • C and ended when temperatures dropped back down to 25 • C. Using this information, for our study we defined the beginning of the summer stressful season as being triggered by three consecutive days with mean water temperatures above 23 • C, and the end of the season being triggered by three consecutive days with mean water temperatures below 25 • C.
Temperature and salinity measurements were taken from a CBNERR continuous water quality station established in 1998 at the Goodwin Islands reserve site. This station houses a YSI (Xylem Analytics) multi-parameter sonde that is fixed to a piling 0.5 m above the bottom that samples every 15-minutes. Data are available through the Virginia Estuarine and Coastal Observing System (VECOS) website [31]. Triplicate water samples were also taken twice a month at this same location as part of a long-term, shallow water monitoring project. Samples were filtered the same day upon return to the laboratory, and then analyzed for chlorophyll [32] and total suspended solids (TSS). Downwelling light attenuation coefficients (K d ) were measured by diffuse downwelling attenuation of photosynthetically available radiation (PAR) determined using triplicate water column measurements of downwelling photosynthetic photon flux density measured with a LI-COR underwater quantum sensor (190SA, LI-Cor, Inc., Lincoln, NE, USA) taken 10 cm below the surface and just above the bottom. Mean monthly water level data were obtained from the U.S. NOAA, National Ocean Survey tide gauge at the Yorktown USCG Training Center.

Data Analyses
All data were analyzed and figures were produced using RStudio v 3.5.0 [33]. To better understand the trend in the long-term record of summer season lengths, a linear model was fit to summer season length as a function of year. To analyze the relationship between seagrass coverage and water temperature, a linear model was fit to seagrass coverage as a function of the length of the summer season. Because seagrass imagery is obtained in June, at the beginning of the summer season, the previous year's length of season was used in the regression.
Water quality parameters were separated into declining and increasing seagrass periods to match the annual seagrass coverage data. Since seagrass imagery is typically obtained in June, we set the cut-off for the periods to start on June 1st and end on May 31st just before the start of the imagery for the next period. For example, for the declining years of 2001-2006, water quality data were used from 01-06-2001 to 31-05-2006. We were interested in differences between the two periods across different months. ANOVA assumptions were not met for equal variances, so Kruskal-Wallis tests were conducted for each month, using the seagrass period (decline or increase) as the factor for salinity, TSS, chlorophyll, K d , surface PAR, and tidal height. A frequency curve for summer (June-August) water temperature was used to compare the two periods, binning the data into the percent time spent at each temperature bin.

Results
Averaged across all years (2004-2017), the maximum biomass of Z. marina follows a unimodal pattern, peaking in May and June and declining by July, with little to no recovery in the fall ( Figure 2). Historical data from this same region from the 1980s and 1990s show Z. marina biomass peaking in July, with declines beginning in August, and some recovery in October. (Figure 2) [17]. Z. marina's historical biomass is more than double that of the more recent surveys, with maximum values of 220 g DW m −2 from surveys conducted in the 1980s and 1990s, compared with a maximum of 102 g DW m −2 from 2004-2017. Investigations of temperature patterns that may be driving these changes revealed a significant increase in the length of the stressful summer season ( Figure 3). Since 1950, the length of this season has been extended by approximately 22 days. This is driven by both an earlier start to the season as well as a later end to the season (data not shown).  Investigations of temperature patterns that may be driving these changes revealed a significant increase in the length of the stressful summer season ( Figure 3). Since 1950, the length of this season has been extended by approximately 22 days. This is driven by both an earlier start to the season as well as a later end to the season (data not shown). Investigations of temperature patterns that may be driving these changes revealed a significant increase in the length of the stressful summer season ( Figure 3). Since 1950, the length of this season has been extended by approximately 22 days. This is driven by both an earlier start to the season as well as a later end to the season (data not shown).   Aerial photography from the Mobjack Bay segment, which includes Goodwin Islands, revealed a significant negative relationship between seagrass coverage and the length of the stressful summer season ( Figure 5). From 1998-2017, there were two distinct periods of time where seagrass hectares (ha) were declining, followed by periods of recovery ( Figure 6). During the first period of decline, seagrasses declined from 3850 ha in 2001 to 2183 ha by 2006, a 43% loss. Recovery occurred, peaking in 2009 at 3521 ha, followed by another decline, reaching a low of 2133 ha in 2012, a loss of 39%. The beds were once again able to recover, peaking in 2015 at 3319 ha.  During multi-year periods of seagrass decline, a greater portion of the time was spent at water temperatures above 28 °C compared with periods of increase ( Figure 7). During the declining periods, the seagrass spent 31% of the time above this threshold, while periods of seagrass increases spent 23% of the time above 28 °C. Mean salinity was typically lower during declining years, ranging between a minimum of 17.2 in April to a maximum of 20.6 in August, while increasing years ranged between a minimum of 17.5 in May to a maximum of 21.8 in November (Figure 8).  During multi-year periods of seagrass decline, a greater portion of the time was spent at water temperatures above 28 • C compared with periods of increase ( Figure 7). During the declining periods, the seagrass spent 31% of the time above this threshold, while periods of seagrass increases spent 23% of the time above 28 • C. Mean salinity was typically lower during declining years, ranging between a minimum of 17.2 in April to a maximum of 20.6 in August, while increasing years ranged between a minimum of 17.5 in May to a maximum of 21.8 in November (Figure 8).

Figure 5.
Linear relationship between the length of the previous summer stressful season and seagrass coverage (hectares). Seagrass hectares were calculated from the VIMS annual aerial mapping survey. Only years where data for both seagrass coverage and the length of the season were available are included. y = −22x + 5706; p < 0.01; r 2 = 0.17. During multi-year periods of seagrass decline, a greater portion of the time was spent at water temperatures above 28 °C compared with periods of increase ( Figure 7). During the declining periods, the seagrass spent 31% of the time above this threshold, while periods of seagrass increases spent 23% of the time above 28 °C. Mean salinity was typically lower during declining years, ranging between a minimum of 17.2 in April to a maximum of 20.6 in August, while increasing years ranged between a minimum of 17.5 in May to a maximum of 21.8 in November (Figure 8).  Water quality comparisons during these time periods revealed overall higher TSS, chlorophyll, and Kd during the declining years compared with the recovery years, though this varied some depending on the month of the year (Figure 9). For the water clarity metrics, differences were greatest in April, where declining years had a mean TSS of 15.4 mg L −1 compared with 9.2 for increasing years, a mean chlorophyll of 14.3 µg −1 compared with 6.7, and a mean Kd of 1.1 m −1 compared with 0.80. Although Kd was greater in April during declining years, the absolute light available did not differ, as mean surface PAR readings were not significantly different. Monthly mean tidal height showed no significant trends, and did not differ among declining vs. increasing seagrass years. Water quality comparisons during these time periods revealed overall higher TSS, chlorophyll, and K d during the declining years compared with the recovery years, though this varied some depending on the month of the year (Figure 9). For the water clarity metrics, differences were greatest in April, where declining years had a mean TSS of 15.4 mg L −1 compared with 9.2 for increasing years, a mean chlorophyll of 14.3 µg −1 compared with 6.7, and a mean K d of 1.1 m −1 compared with 0.80. Although K d was greater in April during declining years, the absolute light available did not differ, as mean surface PAR readings were not significantly different. Monthly mean tidal height showed no significant trends, and did not differ among declining vs. increasing seagrass years.  After the 2010 die-off events, there were very few reproductive shoots present in the spring of 2011 ( Figure 10). During the bed's recovery, there was a sharp increase, so that by spring of 2013, a mean of 53% of the shoots were reproductive. This began to decline through 2015, before increasing again after the 2015 die-off event, leveling out at 20% in the spring of 2017.

Discussion
This study provided evidence for a climate-induced shift in the seasonal growing pattern of a marine foundational species, Z. marina. We related this shift to an increase in the length of the stressful summer season, impacting the timing of Z. marina's growth and decline periods. Over the past 20 years, seagrass growing in the York River and Mobjack Bay areas of the Chesapeake Bay have exhibited distinct periods of decline, characterized by higher water temperatures and lower water clarity compared with the periods of seagrass increase and recovery. Signs of seagrass bed resiliency include an increase in Z. marina reproductive effort following declines, as well as an expansion of a second species, R. maritima, into areas previously dominated by Z. marina.

Growth Pattern
A recent analysis conducted by the University of Maryland's Horn Point Environmental Laboratory found that spring is starting earlier and fall is starting later in the Chesapeake Bay, with the agricultural growing season having increased by more than 30 days over the last century. [34]. When we analyzed the York River temperature data, focusing on a specific definition of the stressful summer season for Z. marina, we found similar results, with the stressful season increasing by 22 days

Discussion
This study provided evidence for a climate-induced shift in the seasonal growing pattern of a marine foundational species, Z. marina. We related this shift to an increase in the length of the stressful summer season, impacting the timing of Z. marina's growth and decline periods. Over the past 20 years, seagrass growing in the York River and Mobjack Bay areas of the Chesapeake Bay have exhibited distinct periods of decline, characterized by higher water temperatures and lower water clarity compared with the periods of seagrass increase and recovery. Signs of seagrass bed resiliency include an increase in Z. marina reproductive effort following declines, as well as an expansion of a second species, R. maritima, into areas previously dominated by Z. marina.

Growth Pattern
A recent analysis conducted by the University of Maryland's Horn Point Environmental Laboratory found that spring is starting earlier and fall is starting later in the Chesapeake Bay, with the agricultural growing season having increased by more than 30 days over the last century [34]. When we analyzed the York River temperature data, focusing on a specific definition of the stressful summer season for Z. marina, we found similar results, with the stressful season increasing by 22 days since 1950. These long-term temperature changes are occurring relatively rapidly, and Z. marina has responded by shifting its growth to peak earlier, a phenomenon documented globally in terrestrial plants [35][36][37]. While a longer growing season may benefit many species through enhanced vegetative growth [38,39], species such as Z. marina that are sensitive to elevated temperatures, particularly when growing near their distributional limits, will be impacted negatively. Z. marina total non-structural carbohydrate reserves peak during the spring and early summer in this region, and the plants can survive the stressful summer by relying on these stored reserves [40]. With a longer stressful season and the plants peaking earlier, they are having to rely on these stored reserves for longer periods of time. Historically, Z. marina in this region exhibited a bimodal growth pattern, with a second peak in the fall. Earlier depletions of these carbon reserves combined with the stressful growing season extending later may be why we more recently see little to no fall regrowth once temperatures cool down [20,41].
Along with temporal shifts, spatial shifts in Z. marina growth patterns have been documented as well. With elevated temperatures increasing light requirements, their mean depth has been declining, and the beds have been moving closer to shore [12,42]. Under a combination of anthropogenic stressors and a changing climate, Z. marina is experiencing habitat squeeze both temporally and spatially, with a shorter period of growth when temperatures are tolerable, and a shrinking depth zone where it is able to get enough light.

Climate and Water Quality Conditions
Water temperatures during declining seagrass years were above 28 • C for a greater portion of the time compared with increasing seagrass years. Seagrass in this area are stressed, with growth rates declining throughout the summer, even at temperatures below this, however this appears to be a critical threshold between plants simply being stressed, and meadows experiencing widespread mortality [21,41].
Clear differences emerged in water clarity metrics when comparing periods of seagrass decline and periods of increase, with the spring and fall periods particularly standing out. April stood out as a month with some of the largest differences among TSS, chlorophyll, and K d . These water clarity metrics appear to dominate, regardless of absolute light availability or tidal height, which did not differ. The cooler water temperatures compared with the summer, and more available light during the longer days compared with the winter, makes this spring period a critical time for Z. marina to accumulate carbon reserves to survive the stressful summer period when respiration exceeds photosynthesis [40].
Interestingly, all light metrics during this critical spring period met or exceeded the Chesapeake Bay habitat requirements [14,15], even for the periods where seagrass was declining. K d values were well within the suitable habitat range across all months, even though seagrass was declining during those periods. These habitat requirements were put in place before climate change and elevated water temperatures were viewed as critical components to Z. marina's growth. This study provides evidence that at least in some instances, these specific water clarity metrics may no longer be applicable in today's climate. With warmer temperatures, Z. marina requires more light to maintain a positive carbon balance [20,24], so with water temperatures continuing to increase, and the length of the stressful summer season increasing, Z. marina will require more light, especially during the critical spring season, to persist. Additionally, following warming events in which significant seagrass habitat is lost, the bed loses its ability to trap sediment, causing greater resuspension and reduced water clarity, creating a negative feedback loop [21]. Salinity was typically lower during declining seagrass years, with the lowest being in April. This corresponds with the relatively lower water clarity also measured during the spring period, indicating seagrass declines are also associated with higher levels of precipitation. Precipitation is projected to increase in the winter and spring for this area, along with streamflow [43]. With water clarity being particularly important during the spring period, this will likely continue to be a major issue for Z. marina.

Resiliency
Seagrasses can exhibit resiliency through four dominant biological features, including genetic diversity, species diversity, continuous habitat, and species biological traits [44]. One of the species biological traits for Z. marina that has played a critical role in its recovery is the reliance on sexual reproduction [19,45,46], existing in an annual form in highly stressful environments [47]. In our study area, Z. marina was able to recover from a heat-induced die-off event in 2005 by a large seed germination event, in which the majority of regrowth the following spring occurred via seedlings rather than clonal growth [19]. Our reproductive shoot data provide evidence for the same strategy following the 2010 die-off event, in which a very small percentage of the shoots in the spring of 2011 were reproductive, presumably due to the fact that the majority of the regrowth was from seedlings that do not produce reproductive shoots in their first year of growth [48,49]. Over the next few years, as the bed was in recovery, the reproductive output was increasing, with reproductive shoot densities reaching 53% of total shoot densities three years after the die-off event. This is higher than reported values for this species in the Chesapeake Bay [48,50] and North Carolina [51]. This increase in reproductive output during the first few years following a stressful event may be an adaptive strategy for this species, spending more resources on seed bank development to ensure survival. Four years into the bed's recovery, reproductive output began to decline to levels more typically seen in this area (11%).
The seagrass beds studied here also exhibited another mechanism for resiliency through species diversity. The interaction between the two species in our system, R. maritima and Z. marina, can be seen as an example of "response diversity", defined as the "diversity of responses to environmental change among species that contribute to the same ecosystem function" [52]. R. maritima's ability to rapidly expand and temporarily inhabit areas where Z. marina has declined appears to be critical for the recovery and continuing functionality of these seagrass meadows. A similar example was documented in San Diego, California following an El Niño event, in which R. maritima temporarily replaced Z. marina after it declined due to elevated temperatures [53]. This temporary replacement may enhance the resiliency of these meadows through positive feedback related to sediment trapping and improvements in water clarity, though the extent to which R. maritima can fill this role is not fully understood. Climate change is likely to continue to favor R. maritima expansion, possibly to the detriment of Z. marina [54]. How these two species interact, and what their interaction means in terms of seagrass bed productivity, will continue to be an important area of research.