Growth and Photoacclimation Strategies of Three Zostera Species Along a Vertical Gradient: Implications for Seagrass Zonation Patterns

Seagrasses typically display vertical zonation patterns in tropical and temperate areas. In East Asia, three Zostera species (Z. japonica, Z. marina, and Z. caespitosa) are distributed at different water depths, leading to zonation patterns from intertidal to subtidal zones. The present study aimed to determine whether these patterns could be explained by differences in their physiological responses, morphological traits, and meadows architecture caused by depth variations. Intertidal Zostera japonica had the lowest shoot height and blade width and the highest shoot density, while subtidal Z. caespitosa exhibited the largest specific leaf area. Total carotenoid content and the total chlorophyll/total carotenoid ratio were significantly higher in Z. japonica than in Z. marina or Z. caespitosa. The average carbon stable isotope (δ13C) values were significantly less negative in Z. caespitosa than in Z. marina and Z. japonica. The average maximum relative electron transport rate was higher for Z. marina than for Z. japonica and Z. caespitosa, while the maximum quantum yield values were higher for Z caespitosa than for Z. japonica and Z. marina. The non-photochemical quenching (NPQ) of Z. japonica extinguished within 5-min following the rapid light curve, indicating that most of its NPQ was energy-dependent quenching (qE). The rapid light curve (RLC) parameters and leaf relaxation after RLC showed that Z. japonica is a high-light adapted species and Z. caespitosa is a low-light adapted one. Additionally, Z. caespitosa showed the lowest dark respiration values and Z. japonica, the highest, suggesting that Z. caespitosa can survive with low photosynthetic activity. These results suggested that interspecific variations in the morphological and physiological characteristics, photosynthetic pigments, and meadows architecture of three Zostera species were associated with vertical zonation patterns and likely explain their patterns.


INTRODUCTION
Seagrasses are a specialized group of flowering plants adapted to marine environments and play crucial functional and structural roles in estuarine and coastal ecosystems (Kuo and McComb, 1989;Orth et al., 2006). They provide food and shelter that support the productivity and biomasses of commercially important fish species (Hemminga and Duarte, 2000;Beck et al., 2001;Short et al., 2007). Additionally, seagrass meadows are noted for their carbon-sink capacity (Nellemann and Corcoran, 2009). However, estuarine and coastal marine environments are undergoing unexpected and drastic changes such as nutrient enrichment, habitat destruction, rising temperatures, light reduction, and sea level rises due to anthropogenic and natural causes (Brierley and Kingsford, 2009;Short et al., 2016). Intertidal and shallow subtidal areas are highly susceptible to such environmental changes; thus, the distribution and abundance of seagrasses in these areas may be particularly vulnerable to such changes (Orth et al., 2006). Furthermore, climate change greatly accelerates environmental changes and, consequently, causes considerable changes in seagrass distributions (Duarte et al., 2018). Thus, understanding how seagrass species respond to environmental changes is critical for developing effective strategies to protect seagrass beds.
Seagrasses usually form single-species meadows or mixed vegetation areas with high species diversity in estuaries and lagoons and show specific zonation patterns from intertidal to subtidal areas in the temperate seagrass meadows (Short et al., 2007). Biotic factors such as grazing and interspecific competition affect seagrass zonation (Tribble, 1981;Williams, 1987). The morphological and physiological responses of seagrasses to local environmental conditions and stressors are well-studied (Zhang et al., 2015;Schubert et al., 2018;Soissons et al., 2018). For example, intertidal seagrasses exposed to excess light, desiccation stresses, and a wide temperature range had narrower leaves and a shorter canopy (Manassa et al., 2017). The increased shoot density of intertidal seagrass may be an adaptation to minimize emergence or desiccation stress (Yabe et al., 1996). Contrastingly, subtidal seagrasses increase the width and length of their leaves to enhance light capture under low light availability (Lee and Dunton, 1997;Collier et al., 2007). Moreover, seagrasses in intertidal and subtidal zones showed clearly different physiological responses to light conditions (Ralph et al., 1998;Park et al., 2016). These zonation patterns reflect their morphological and physiological adaptions to local environmental conditions and stressors, or their ability to occupy various niches (Björk et al., 1999;Tanaka and Nakaoka, 2004;Kaldy et al., 2015).
Intertidal seagrass species may be more exposed to excess irradiance, wide temperature ranges and desiccation stress than subtidal seagrass species due to tidal changes. The distribution limits of intertidal seagrass species may be determined by their adaptation or tolerance to air exposure, interspecific competition, light availability, and physical disturbance from tidal currents (Boese et al., 2005;Cabaço et al., 2009;Van der Heide et al., 2010). Leuschner et al. (1998) highlighted that the difference in desiccation tolerance between Zostera marina and Z. noltii allows them to have different depth distributions. Conversely, Z. japonica distributions in North America can be explained by temperature rather than desiccation tolerance and light limitation (Kaldy et al., 2015). In submerged seagrass beds, light availability determines seagrass depth distributions (Duarte, 1991). Since seagrasses require more light for growth and survival than other marine macrophytes, determining the threshold light levels for seagrass growth and survival is crucial for efficiently managing seagrass beds (Björk et al., 1999). For example, seagrasses respond to light intensity reductions by increasing their chlorophyll content, decreasing their chlorophyll a/b ratio and reducing their production of ultraviolet (UV) lightblocking pigments to improve light capture and chemical energy conversion (Abal et al., 1994;Park et al., 2016). Additionally, the physiological responses of seagrasses to reductions in light intensity vary with temperature and, therefore, season as well (Collier et al., 2008;York et al., 2013). As the factors contributing to vertical zonation patterns show seasonal variations, it is essential to investigate the seasonal responses of seagrasses to evaluate the depth-related differences in their morphology and physiology.
Fifteen seagrass species are distributed in the temperate North Pacific region, where three genera (Zostera, Ruppia, and Phyllospadix) dominate (Short et al., 2007). In particular, the vertical patterns of five Zostera species (Z. japonica, Z. marina, Z. caespitosa, Z. caulescens, and Z. asiatica) in Northeast Asia (China, Japan, and Korea) occur along waterdepth gradients (Aioi and Nakaoka, 2003;Lee and Lee, 2003). While Z. caulescens and Z. asiatica occur in relatively deep zones and form single-species meadows, Z. japonica, occurs in the intertidal zone, Z. marina grows to about 10 m depth, and Z. caespitosa reach 20 m depth in the temperate North Pacific region (Short et al., 2007). Additionally, three Zostera species form mixed vegetation areas (Lee et al., 2018). The intraspecific variations in the morphological and physiological characteristics of seagrasses along depth variations are relatively well-documented (Olivé et al., 2013;Park et al., 2016;Manassa et al., 2017). Contrastingly, the interspecific variations in their morphological and photosynthetic responses to environmental factors with respect to vertical zonation patterns are unclear (Huong et al., 2003;Tanaka and Nakaoka, 2007;Apichanangkool and Prathep, 2014;Shafer and Kaldy, 2014). Thus, further in situ studies and information on seagrass ecology and physiology are required to improve seagrass ecosystem management.
Here, we aimed to determine whether the observed vertical zonation patterns of three Zostera species could be explained by seasonal variations in physiological responses and morphological traits caused by depth variations. We hypothesized that morphological, meadow architectural, and physiological features of these species reflected depth-related environmental stressors such as light availability, temperature, and air exposure. This was done by examining their morphology, density, biomass, productivity, photosynthetic pigments, carbon stable isotope ratio, chlorophyll a fluorescence parameters, and respiration rates. We also assessed light-dark relaxation kinetics to evaluate photosystem recovery upon high light exposure. Global climate change is a key driver of accelerating environmental change and consequently results in considerable changes in seagrass zonation. Thus, understanding how the three Zostera species respond to changing environmental conditions will help to predict their vertical zonation patterns in the temperate North Pacific region.

Study Area
The study site was located in Koje Bay (34 • 48 ′ N, 128 • 35 ′ E) on the southern coast of Korea, where three Zostera species (Z. japonica, Z. marina, and Z. caespitosa) are distributed sequentially from the intertidal zone to the shallow subtidal zone (Figure 1). The seagrass beds in Koje Bay consist of three monotypic meadows (Z. japonica, Z. marina, and Z. caespitosa) and two transition zones (Z. japonica/Z. marina and Z. marina/Z. caespitosa). Zostera japonica meadows were distributed from the upper to lower intertidal zones (75-45 cm above the mean lower low water, MLLW) and were exposed to air at low tide. Zostera marina meadows grew in the lower intertidal (the lower distributional limit of Z. japonica) and upper subtidal (the upper distributional limit of Z. caespitosa) zones. In contrast, the Z. caespitosa beds occurred at an average water depth of ∼3.0-6.0 m relative to the MLLW. This study was performed in monotypic meadows of Z. japonica, Z. marina, and Z. caespitosa. Zostera japonica and Z. marina meadows were located at 63-67 cm and 15-20 cm above the MLLW, respectively, while Z. caespitosa were located at average water depths of ∼4.5-5.0 m relative to the MLLW. Water temperatures in Koje Bay ranged from 5.0 to 28.5 • C, and salinity ranged from 30.4 to 33.6 (Kim et al., 2015(Kim et al., , 2016aSuonan et al., 2017). Sandy sediments predominated at the study site. The tidal regime in the region is mixed semi-diurnal with a 3.0 m tidal range (Tide tables for the Coast of Korea, National Oceanographic Research Institute of Korea). The most extreme tidal range occurs in early spring, while the smallest tidal range in Summer. The summer monsoon rainy season and typhoons occur from mid-June through late-September.

Light Irradiance and Exposure Time
Light irradiance was monitored every 15 min using HOBO light intensity data loggers (Onset Computer Corp., Bourne, MA, USA) enclosed in clear submersible polycarbonate cases at the three study zones. The underwater sensors were cleaned regularly to minimize fouling. The light intensity (lumens ft −2 ) measured using the HOBO data loggers was converted to photon flux density (PFD; mol photons m −2 s −1 ) by concurrent quantum measurements using an LI-1400 data logger and an LI-193SA spherical quantum sensor (LI-COR, Lincoln, NE, USA). The HOBO data were calibrated to the LI-COR data using an exponential decay fit (Long et al., 2012). The HOBO measurements and PFD exhibited a relatively strong correlation (r 2 = 0.90). Daily PFD (mol m −2 d −1 ) was calculated as the sum of quantum flux over each 24-h period.
The tidal height at each zone was measured when the study zone began to be exposed to the air during low tide, and the air exposure time was calculated at each station using real-time tidal height data for the Korean coastline (Korea Hydrographic and Oceanographic Administration; http://www.khoa.go.kr).
Monthly exposure time (h) was represented as the sum of total exposure time in each month.

Morphology and Pigments
To measure shoot morphological characteristics, 10-15 mature terminal shoots from each species were collected monthly from November 2005 to September 2006. Shoot height, sheath length and the width of the longest leaf were measured to the nearest 0.1 mm. To determine the specific leaf area (SLA; cm 2 g −1 dry weight), leaf samples of 1 cm 2 were placed in a drying oven at 60 • C to a constant weight. All leaves and below-ground tissues of the shoots were dried at 60 • C to a constant weight to determine their individual shoot weight.
To determine blade chlorophyll content, six replicate samples from each bed were collected and cleaned in the laboratory. Leaf tissues were extracted at 4 • C for 2-3 days in glass tubes containing 5.0 mL of N,N-dimethylformamide (DMF) in the dark before spectrophotometric determination (Shimadzu UV model 2264). Absorbance was measured at 480, 647, and 664 nm and the absorbance at 750 nm was subtracted to correct for scattering due to turbidity. Pigment contents were calculated using the Wellburn (1994) extinction coefficient equations.

Shoot Density, Biomass, and Production
Shoot density inside randomly thrown quadrats or cores (Z. japonica, 10-cm diameter core; Z. marina and Z. caespitosa, 0.35 × 0.35 m; n = 4-6) were counted bimonthly to minimize the damage to seagrass beds and were converted to per unit area (shoots m −2 ). Shoot biomass was calculated by multiplying shoot density by the average individual weight per shoot and was converted to per-unit-area estimates (g dry weight [DW] m −2 ). Leaf productivity was measured using a modified blademarking technique (Zieman, 1974;Kentula and McIntire, 1986). Ten to fifteen randomly chosen shoots were marked through the sheath bundle above the meristem using a hypodermic needle and they were harvested after 2-4 weeks (Park et al., 2010). Leaf materials were separated into the leaf tissues produced before and after marking, and they were then dried at 60 • C to a constant weight. Leaf productivity per shoot (mg DW shoot −1 d −1 ) was determined by dividing the dry weight of the new leaf tissues produced after marking by the number of days since marking. Areal leaf productivity (g DW m −2 d −1 ) was calculated by multiplying productivity per shoot by shoot density.

Stable Isotope Analysis
The youngest mature leaf (second or third leaf) of the shoots collected for the morphological analysis was used to assess the carbon stable isotope ratio (δ 13 C) of leaf tissues. Leaf tissues were dried at 60 • C to a constant weight and were ground using a mortar and pestle. The ground samples were wrapped in tin capsules and were analyzed to determine their carbon stable isotope ratios using an elemental analyzer (Eurovector 3000 Series; Milan, Italy) coupled with a continuous-flow isotope ratio mass spectrometer (Isoprime, GV Instruments). For carbon, Pee Dee Belemnite marine limestone (PDB) was used as the primary standard. Analytical precision was ∼0.1‰ for δ 13 C.

Chlorophyll a Fluorescence Measurements
The chlorophyll a fluorescence of the three seagrass species (n = 10-15) were measured in situ using a pulse amplitude modulated (PAM) fluorometer (Diving-PAM, Walz) with an 8-mm optical fiber. The sampling position of the optical fiber was at the midpoint of the youngest mature leaf and the tip of the instrument's optical fiber was placed 6 mm from, and perpendicular to, the middle of the adaxial surface of the leaves using dark leaf clips (DLC-8). The maximum quantum yield (F v /F m ) of photosystem II (PSII) was determined after 10-min of dark-adaptation using leaf tips and was calculated using the following equation: where F 0 is the minimal fluorescence of a darkacclimated leaf in which all PSII reaction centers are open, and F m is the corresponding maximum fluorescence measured with all PSII reaction centers closed following a short period (0.8 s) of saturating light (e.g., Van Kooten and Snel, 1990).
The rapid-light curve (RLC) was produced by the Diving-PAM fluorometer under an internal program using artificial photosynthetic photon flux density (PPFD). The samples were exposed to each irradiance level (56,115,184,268,366,537,721, and 1042 µmol photons m −2 s −1 ) at 10-s intervals. The RLCs used a measure of relative electron transport rates (rETRs) determined from the following equation: rETR = F/F m ' × PPFD × 0.5 × 0.84, where F/F m ' is the effective quantum yield of PSII, 0.84 is the assumed absorption coefficient, 0.5 is the correction for two photosystems absorbing photons, and PPFD is the programmed level of photosynthetically active radiation (PAR, 400-700 nm) delivered by the halogen lamp (Genty et al., 1989;Ralph and Gademann, 2005). A standard absorption factor (0.84) was used because the best means of directly measuring absorption continues to be debated among the scientific community (Collier et al., 2008). As such, all values were considered relative ETRs (rETRs) because leaf absorptance was not directly measured (Durako et al., 2003). The RLC data were fitted to the model of Platt et al. (1981) to obtain values for photosynthetic efficiency (α), inhibition term (β), minimum saturating irradiance (E k ), and maximum relative electron transport rate (rETR max ). E k was calculated as E k = rETR max /α.

Light-Dark Relaxation Kinetics
To assess leaf relaxation after light exposure, the recovery of F/F m ' and other fluorescence parameters that were in darkness for 10-min after an RLC were measured. The RLC + recovery (LC+REC) was measured in the dark at 10 s, 30 s, 1-min, 2-min, 5-min, and 10-min after the last actinic light period of an RLC. The effective quantum yield ( F/F m ') and nonphotochemical quenching (NPQ) parameters were calculated using the following equations:

Respiration Measurements
Dark respiration of entire plants was measured using the Fibox 3 spot sensor and oxygen dipping probe (PreSens, Regensburg, Germany). Seagrass samples with rhizome and root tissues were collected in March 2008 for respiration measurements and were washed in seawater to remove sediment and epiphytes.
Before determining the effects of temperature on respiration, the plants were allowed to equilibrate overnight to the experimental temperature. All experiments were conducted within 2 days of collection. The plants were placed in 500 mL of BOD (biological oxygen demand) bottles filled with 0.2-µm pore sized filtered seawater equipped with magnetic stirrers. The bottles were incubated in the dark in treatment chambers at 8, 16, 22, and 27 • C (average seasonal water temperatures) for 3 h. Next, the plants were removed from the bottles, which were stirred for 15-min to equalize the oxygen concentration in the bottles. The oxygen concentrations were measured before and after incubation for 3 h. In the control incubations (seawater only), there were almost no changes in oxygen concentrations over the 3 h period. Respiration rates were normalized to the dry weight of the plant tissues in the bottles and were expressed as µmol O 2 g −1 DW h −1 .

Statistical Analyses
All values were expressed as means ± standard error (SE). Statistical analyses were performed using SPSS ver. 18.0 (SPSS Inc., Chicago, IL). Data were tested for normality and homogeneity of variance to meet parametric test assumptions. If these assumptions were not satisfied, the data were log or squareroot transformed. Differences in light irradiance, morphological characteristics, photosynthetic pigments, shoot densities, biomasses, leaf productivity, photosynthetic characteristics, and δ 13 C values among the sampling times and species were tested for significance using a two-way analysis of variance (ANOVA), with time as a block. Dark respiration was also analyzed using a two-way ANOVA (species × water temperature). The Student-Newman-Keuls (SNK) test was used to identify statistically significant differences. Statistical significance was set to alpha < 0.05.

Light Irradiance and Exposure Time
Average daily light availability showed clear seasonal variations and was significantly (P < 0.001) different among three zones; it was higher during boreal winter and spring than during summer and fall (Figure 2A). Average daily light irradiances at the three study zones were ∼40.6 (Z. japonica zone), 26.5 (Z. marina zone), and 14.7 mol photons m −2 d −1 (Z. caespitosa zone), respectively (Figure 2A). Monthly air exposure time varied with three zones, being highest during winter-spring and being lowest during summer-fall due to seasonal tidal range, rainy season and typhoons. Zostera caespitosa meadows were never exposed to air, while Z. japonica and Z. marina were exposed to air for approximately bout 2.5 and 0.5 h per day, respectively (Figure 2B).

Morphology and Pigments
Shoot height, sheath length, leaf width, and specific leaf area (SLA) of the three Zostera species exhibited clear and differing seasonal variations (Figure 3). Shoot height, sheath length, and leaf width of Z. marina and Z. caespitosa increased during spring and decreased during fall, while those of Z. japonica were lowest during February-April, when the duration of air exposure was longest. The SLA of the three species increased during winterspring and decreased during late spring-summer. Shoot height, sheath length, and leaf width were highest in Z. caespitosa and lowest in Z. japonica; however, Z. japonica showed the highest SLA and Z. caespitosa, the lowest. Zostera marina exhibited shoot height, sheath length, leaf width, and SLA values that ranged between those of the other two species.
Total chlorophyll (chl. a+b) and total carotenoids contents were highest during summer-fall and lowest during springsummer (Figure 4). Zostera japonica had significantly (P < 0.001 for both) higher total chlorophyll and total carotenoids contents than Z. marina and Z. caespitosa. The chlorophyll/carotenoids ratio varied with sampling time but did not exhibit a clear seasonal trend. The average chlorophyll/carotenoids ratio was lowest in Z. caespitosa and highest in Z. japonica and Z. marina.
Shoot Density, Biomass, Productivity, and δ 13 C Shoot density of the three species exhibited clear seasonal variations, increasing during spring and decreasing during fall and winter (Figure 5). Shoot density of Z. japonica (4,498 shoots m −2 ) was significantly (P < 0.001) higher than that of Z. marina (1,144 shoots m −2 ) and Z. caespitosa (95 shoots m −2 ). The total biomass of the species showed seasonal trends similar to those of shoot density and was significantly (P < 0.001) higher in Z. marina (393 ± 55 g DW m −2 ) than in Z. japonica (123 ± 26 g DW m −2 ) and Z. caespitosa (121 ± 21 g DW m −2 ) ( Figure 5).
Leaf productivity per shoot and areal leaf productivity also exhibited clear seasonal variations, increasing during spring and decreasing during summer and early fall (Figure 6). Zostera caespitosa showed significantly (P < 0.001) higher leaf productivity per shoot than Z. japonica and Z. marina. In contrast, Z. marina showed the highest in areal leaf productivity. The average areal leaf productivity values of Z. marina, Z. japonica, and Z. caespitosa were 5.4 ± 2.1, 1.4 ± 0.4, and 1.0 ± 0.4 g DW m −2 d −1 , respectively.
The δ 13 C values significantly (P < 0.001) varied with species and were dependent on their depth (Figure 7). The δ 13 C values of both Z. japonica and Z. marina peaked in March and decreased in late spring, while those of Z. caespitosa were relatively constant during the experimental period, except in July. The average δ 13 C values were significantly (P < 0.001) more negative in Z. japonica (−10.41 ± 0.42‰) than in Z. marina (−9.33 ± 0.34‰) and Z. caespitosa (−7.57 ± 0.27‰).

Fluorescence Parameters and Dark Respiration
Significant seasonal differences (P < 0.001) in rETR max were observed for Z. japonica and Z. marina, which were highest in January and lowest in July and September, while those of Z. caespitosa were relatively constant ( Figure 8A). The average rETR max values of Z. marina, Z. japonica, and Z. caespitosa were 17.1 ± 1.6, 13.8 ± 1.1, and 11.5 ± 0.4 µmol electrons m −2 s −1 , respectively. The F v /F m and photosynthetic efficiency (α) of the three species exhibited distinct seasonal variations, increasing during winter and spring and decreasing during summer and fall (Figures 8B,C). The F v /F m was highest in Z. caespitosa and lowest in Z. japonica. The photosynthetic efficiency (α) of Z. japonica was significantly (P < 0.001) lower than that of the other two species. The minimum saturating irradiance (E k ) also showed significant (P < 0.001) seasonal variations, but with the inverse trend ( Figure 8D). Zostera japonica showed the highest E k and was significantly different (P < 0.001) from that of the other two species.
During the course of the RLC, the NPQ of Z. japonica was two-to three-fold higher than that of Z. marina and Z. caespitosa (Figure 9). The relative relaxation proportion of Z. japonica (56%) was lower than that of Z. marina (78%) and Z. caespitosa (71%) during the first 30 s. However, the relative relaxation proportion of Z. japonica (98%) was significantly (P < 0.001) lower than that of Z. marina (93%) and Z. caespitosa (87%) after 2-min. The NPQ of Z. japonica was completely extinguished after 5-min (100%) in the dark. In contrast, even after 10-min of darkness, around 1.4 and 5% of NPQ remained in Z. marina and Z. caespitosa leaves, respectively.
Dark respiration in the three species varied significantly (P < 0.001) with water temperature ( Table 1). The average dark respiration of Z. japonica was significantly (P < 0.001) higher than that of Z. marina and Z. caespitosa (Table 1). Additionally, the dark respiration of the three species increased significantly with increasing water temperature (P < 0.001 in species × water temperature interaction), which was highest at 27 • C and lowest at 8 • C.

DISCUSSION
We found interspecific variations in the morphological and photosynthetic characteristics of three Zostera species associated with vertical zonation patterns. Subtidal Z. caespitosa was the tallest of the three species and had wider leaves, while Z. japonica in the upper intertidal zone was the shortest and had the thickest and narrowest leaves. The morphological responses of seagrasses to reduced light availability with increasing depth have been widely reported (Olesen et al., 2002;Collier et al., 2007;Mackey et al., 2007;Kohlmeier et al., 2014). In shading experiments with Z. capricorni, shoots exposed to high light irradiance were shorter than those exposed to low light irradiance (Abal et al., 1994). Subtidal seagrasses enlarge their shoots or produce longer and wider leaves to capture photons under lowlight conditions (Watanabe et al., 2005;de los Santos et al., 2010). Enríquez and Pantoja-Reyes (2005) reported that leaf morphology (leaf thickness or specific leaf area) is related to the depth-induced reductions in surface irradiance. This shows that light availability influences morphological variations (leaf length and width) in seagrasses and may be the main factors regulating the vertical zonation of seagrasses. Seagrasses have the ability to photo-acclimate by altering and optimizing photosynthetic pigment content in response to light availability (Enríquez, 2005;Ralph et al., 2007). With increasing depth or decreasing light irradiance, leaf chlorophyll a+b content usually increases, while UV absorbance and total carotenoids content typically decreases (Longstaff and Dennison, 1999;Olesen et al., 2002;Park et al., 2016). However, the opposite effect occurs under low-light conditions (Biber et al., 2009;Collier et al., 2012;Kohlmeier et al., 2014). In the present study, chlorophyll and total carotenoids contents, and the chlorophyll/carotenoids ratio of three Zostera species corresponded to ambient light availability. Zostera japonica exposed to more light during low tide had the highest carotenoids content and subtidal Z. caespitosa had a relatively low carotenoids content throughout. Carotenoids lose  excitation energy through heat, which protects photosystems or chlorophyll from excess light (Amstrong and Hearst, 1996;Demmig-Adams and Adams, 1996) and intertidal Z. marina populations have higher carotenoids contents than subtidal populations to avoid photodamage due to excess light during the day (Park et al., 2016). When plants are exposed to high light irradiances, their total carotenoids content increases (Härtel and Grimm, 1998;Ralph et al., 2002). Therefore, the observed differences in pigment content appear to be an important acclimation strategy across all three species enabling them to cope with gradients in light conditions occurring along the depth zonation.
Seagrass vertical zonation patterns can be explained by mechanisms at the population and meadow level and by physiological characteristics at the individual level (Tanaka and Nakaoka, 2004;Shafer et al., 2007). Subtropical seagrasses mitigate desiccation stress by allowing their leaves to lie atop moist sediment (Björk et al., 1999;Tanaka and Nakaoka, 2004) and the high shoot density of intertidal seagrass can reduce desiccation stress through shading and overlapping leaves (Tanaka and Nakaoka, 2004;Shafer et al., 2007). Although the effective quantum yield of Z. japonica at the individual level was more affected by desiccation than that of Z. marina, the effective quantum yield of the densely overlapping leaves of Z. japonica was significantly higher than that of Z. marina with air exposure (Kim et al., 2020). In this study, Z. japonica may have been weakly affected by desiccation stress during exposure at low tide, possibly because its structural flexibility and high shoot density allows leaves to contact moist sediment and reduce water loss. This indicates that growth strategies should be considered when evaluating seagrass desiccation tolerance.
The rETR max and F v /F m values from the RLCs indicated physiological differences and acclimation along the depth-related light gradient. Subtidal Z. caespitosa had low photosynthetic activity, while Z. japonica and Z. marina exposed to air at low tide had high photosynthetic activity. However, the rETR max of Z. japonica was lower than that of Z. marina even though the exposure time of Z. japonica was two-fold that of Z. marina. Prolonged exposure to high irradiance may reduce F v /F m in marine plants (Yakovleva and Titlyanov, 2001;Durako and Kunzelman, 2002;Manassa et al., 2017). Zostera caespitosa that were not exposed to air showed a high F v /F m value, which is consistent with past reports of increase in leaf fluorescence ratios in response to decreased light irradiance in seagrass (Dawson and Dennison, 1996;Major and Dunton, 2002). The photosynthetic efficiency (α) and E k values of the three species also showed consistent patterns with depth. Zostera japonica had low α and high E k values, while Z. caespitosa showed the reverse trend. Zostera marina had intermediate α and E k values. This was similar to seagrass zonation patterns, suggesting successful adaption to environmental conditions. The NPQ can be divided into three processes based on the time required to relax after exposure to high-light conditions. The NPQ components that relax within a range of seconds to minutes are associated with the removal of energy-dependent NPQ (qE) and state-transition NPQ (qT), whereas slower relaxation (>10min to several hours) is associated with photoinhibition (qI) and changes in energy distribution in favor of PSII (Horton and Hague, 1988;White and Critchley, 1999). LC+REC helps distinguish the relaxation kinetics of the NPQ components (Ralph and Gademann, 2005). Although the NPQ of Z. caespitosa was mostly relaxed, ∼5% of the NPQ remained in Z. caespitosa leaves after 10-min of darkness. This is likely to be photoinhibitory quenching (qI). However, the NPQ of Z. japonica was energy-dependent quenching (qE), as it was mostly relaxed after 2-min and completely extinguished after 5-min. This results was similar to that of Ralph and Gademann (2005), who reported that The NPQ in high-light adapted seagrass leaf completely extinguished after 5-min. This suggests that Z. japonica may be a high-light adapted species, while Z. caespitosa may be a low-light adapted species. Shallow-water seagrasses have a higher capacity for light protection and NPQ compared to deep-water seagrasses (Ralph et al., 1998).
Seagrass respiration is primarily influenced by temperature and depth (reviewed by Touchette and Burkholder, 2000). Respiration rates continue to increase with increasing temperatures (Bulthuis, 1983;Marsh et al., 1986;Pedersen et al., 2016). Dark respiration of Z. marina and Thalassia testudinum was lowest in winter, increased with increasing water temperature in spring and peaked in summer (reviewed by Lee et al., 2007). The respiration rate of Z. marina at 5.5 m was 40% lower than that at 1.3 m (Dennison and Alberte, 1982). To prevent anaerobic respiration in below-ground tissues, many aquatic plants (including seagrasses) translocate and release O 2 in the rhizosphere during active photosynthesis (Sand-Jensen et al., 1982;Smith et al., 1984;Caffrey and Kemp, 1991). Therefore, seagrasses with high below-ground biomasses have high respiratory demand (Tanaka and Nakaoka, 2007). In this study, dark respiration increased with increasing water temperature, and subtidal Z. caespitosa showed the lowest respiration rate. Additionally, Z. caespitosa had the lowest below-/above-ground biomass (data not shown), which may have affected its respiration rate. These results suggest that Z. caespitosa has the potential to grow in the subtidal zone despite its low photosynthetic activity.
Seagrass δ 13 C values significantly correlate with seagrass production rates (Cooper and DeNiro, 1989;Hemminga and Mateo, 1996). Most seagrasses use HCO − 3 for photosynthetic carbon reduction, which is the most abundant form of dissolved inorganic carbon in seawater and has a less negative δ 13 C value than dissolved CO 2 due to a low supply of CO 2 in seawater (Hemminga and Mateo, 1996;Beer and Rehnberg, 1997;Invers et al., 2001). Therefore, seagrass δ 13 C values increase during high productivity periods by utilizing HCO − 3 with a less negative δ 13 C value than dissolved CO 2 due to reduced carbon discrimination (Kim et al., 2014). Additionally, the δ 13 C values become more negative with less light due to a reduced photosynthetic and production demand for carbon, allowing more carbon discrimination (Cooper and DeNiro, 1989;  Abal et al., 1994;Hemminga and Mateo, 1996). Contrastingly, the reduced photosynthetic discrimination against the heavier isotope may result in increased δ 13 C values due to photosynthetic carbon limitation during periods of high light (Fourqurean et al., 2005;Campbell and Fourqurean, 2009). In this study, however, an increase in δ 13 C values was not observed during high productivity periods. The δ 13 C values were highest in Z. caespitosa and lowest in Z. japonica, regardless of light irradiance. The productivity-δ 13 C relationship was not observed in the shallow seagrass meadows as carbon isotopic discrimination in seagrass saturates earlier than productivity with increasing light irradiance (Mateo et al., 2001). During the experiment, the δ 13 C values were not controlled by the light irradiance in each zone, which suggests that δ 13 C values are not determined by productivity and light conditions.
Seagrass δ 13 C values are mainly determined by inorganic carbon sources (i.e., HCO − 3 or dissolved CO 2 ) and availability during photosynthesis (Raven et al., 2002;Lepoint et al., 2003). They are also controlled by atmospheric CO 2 use and inorganic carbon concentrating mechanisms, such as CO 2 diffusion, active CO 2 transport and active HCO − 3 influx under different environmental conditions (Raven et al., 2002). The heavier δ 13 C values of tropical seagrass are attributed to the use of HCO − 3 due to a low supply of CO 2 in seawater (Hemminga and Mateo, 1996). The lighter δ 13 C of Halophila ovalis leaves suggests that CO 2 was so abundant that H. ovalis did not take up HCO − 3 in the seawater (Yamamuro et al., 2001). Our results showed that the δ 13 C values of the species were consistent with their zonation patterns. For photosynthesis, subtidal Z. caespitosa may have absorbed inorganic carbon from the water column via active CO 2 transport and HCO − 3 influx, resulting in higher δ 13 C values. However, the lighter δ 13 C values of intertidal Z. japonica may have resulted from atmospheric CO 2 use or a sufficient CO 2 supply in seawater due to air mixing during tidal changes (Kim et al., 2016b). Contrastingly, the values of Z. marina may have been influenced by both processes. Thus, these results suggest that three Zostera species utilize different inorganic carbon sources and carbon concentrating mechanisms for photosynthesis to co-exist or avoid resource competition.

CONCLUSION
Although environmental conditions were perhaps different in this study area, these conclusions are still supported for data collected 15 years ago and providing support for the mismatch between the respiration data measured in seagrasses 3 years after other measurements. In this study, these results strongly support our hypothesis, namely, that morphological and physiological characteristics and meadows architecture of three Zostera species strongly reflected vertical zonation patternsrelated environmental stressors (Figure 10). The photosynthetic activities suggest that Z. japonica is a high-light adapted species and Z. caespitosa is a low-light adapted one, while Z. marina shows high plasticity of growth and physiology to allow this species to occur in the intertidal and subtidal zones (Figure 10). If global climate change and its effects are accelerating over time, this change could lead to increasing space competition among seagrass. Tidal inundation and reduced light availability caused by sea level rise may facilitate the penetration of subtidal species into middle and lower intertidal zones. In contrast, these factors appear to limit expansion of intertidal seagrass species meadows lower boundary to the middle and lower intertidal zones. The distribution of intertidal seagrass species may be shifted upper in the intertidal zone.
However, seagrasses vary in their sensitivity to changes in environmental conditions such as light reduction and temperature; for example, Halophila species occurs in both highand low-light regions, while Thalassia and Cymodocea species are found only in high-light area (Carruthers et al., 2002). Under such circumstances, Thalassia and Cymodocea species is expected to shift upper intertidal zone to achieve optimum lighting conditions. Thus, a better understanding how seagrass species respond to environmental changes, which are being exacerbated by climate change, will provide valuable information for predicting their vertical zonation patterns and designing effective management strategies for seagrass restoration.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
SP and K-SL analyzed the data and wrote the draft. KM and SK conducted the fieldwork and analyzed the data. All authors generated the hypotheses, designed the study, helped improve the manuscript, and approved the submission.