In-Situ Estimates of Net Ecosystem Metabolisms in the Rocky Habitats of Dokdo Islets in the East Sea of Korea

<jats:p>We measured oxygen (O2) fluxes in two major shallow subtidal benthic habitats (kelp bed (KB) and bare rock (BR) covered with crustose coralline algae) of Dokdo islet in the East Sea by applying noninvasive in-situ aquatic eddy covariance (AEC). The AEC device allows time series measurements (~24 h) of three-dimensional velocity (u, v, and w components) and high-resolution dissolved O2. This allows estimation of O2 exchange flux via benthic habitats. Local flow rates and irradiance levels were found to be major factors controlling O2 exchange flux in the rocky habitats. Gross primary production rates tended to be significantly higher in KB (163 mmol O2 m−2 d−1) than in BR (51 mmol O2 m−2 d−1). The net ecosystem metabolisms were assessed as opposite types, with 8 mmol O2 m−2 d−1 in KB (autotrophy) and –12 mmol O2 m−2 d−1 in BR (heterotrophy). Our results indicate that kelp beds are important for organic carbon cycling in rocky coastal waters and that AEC application to macroalgae habitats is a useful assessment approach.</jats:p>

section_content: Macroalgal-dominated habitats, mainly kelp beds (forest), cover about 25% of the world's coastlines, from temperate to subpolar regions in both hemispheres [1] [2] [3]. Kelp beds are one of most productive ecosystems in coastal waters, where they also provide complex biogenic habitats for juveniles and enhance secondary production [1, 2]. Organic matter produced by benthic communities in high-productivity kelp beds can become a food source for adjacent coastal to deep-sea ecosystems [1]. Furthermore, because of their high capacity for carbon fixation and storage, through their physiological cycle, these ecosystems are recognized as a blue carbon sink (i. e., carbon stored in coastal and marine ecosystems) in the context of global carbon dioxide cycles. As a result, kelp bed communities may offer an important bioengineering opportunity in coastal waters [1]. However, increasing anthropogenic pressures (e. g., eutrophication, overfishing, climate change) have caused a decline of 38% in the kelp populations found in kelp beds [1, 2]. 

Rocky bottoms, in the shallow euphotic zone of Korea, are important habitats for macroalgae: predominantly, long-lived canopy-forming Ecklonia cava and Eisenia bicyclis [4, 5]. The distribution and abundance of these species are rapidly declining in the coastal waters of Korea, similar to Japan [6]. Although the reasons for the destruction of the kelp forest on the northeast Pacific coasts are not well known, warming of surface water temperatures by Kuroshio Current may be an important key to controlling the distribution of kelp [6]. In addition, a massive outbreak of crustose coralline algae (CCA) whitening has been rapidly expanding on the rocky bottom [7, 8]. These issues can represent a severe threat to the macroalgal communities, but impacts on ecological functions in coastal waters are still limited. 

Production and respiration in the benthic community are key parameters for understanding ecological health and organic carbon flow in the food web [9]. Net ecosystem metabolism (NEM), the difference between gross primary production (GPP) and respiration (R) (NEM = GPP -R), is a valuable proxy for several things: (1) consumption through food webs, (2) remineralization/decomposition by microbial respiration, (3) burial to sediment, and (4) exports from vegetated habitats [9]. Therefore, quantification of NEM is an important component, essential for understanding the ecological functioning of kelp beds. 

Several methods for NEM estimation have been used to date. A simple laboratory approach, based on the conventional light/dark incubation method, is widely used to quantify NEM in benthic communities [10]. The recently developed in-situ photorespirometer overcomes the limitations of laboratory methods [11, 12]. However, the installation of a rigid enclosed chamber in the field can isolate seaweed samples from natural water movements, which can alter the diffusive boundary layer at the thallus surface [13, 14]. Thus, an alternative and robust method that can include surrounding physical effects is preferred to measure NEM in benthic habitats. 

Aquatic eddy covariance (AEC) has been developed in recent decades [15]. AEC is a noninvasive in-situ measurement technique that has been used to measure O 2 flux in various natural water habitats (e. g., coral reefs, seagrass meadows, hard bottom benthic flora beds, oyster reefs) [16] [17] [18] [19] [20] [21] [22] [23] [24]. The measurement covers a relatively large surface area (10-50 m 2 ) and interferes minimally with the prevailing hydrodynamic conditions [25]. In particular, AEC can be validly applied in complex hard-bottom areas, such as rocky bottoms, coral reefs, and oyster reefs, where conventional chamber deployment is problematic [17, 18, 21, 23]. 

In this study, we estimated, in situ, the metabolism balance of benthic communities in vegetated (kelp bed [KB]) and 'bare' rock [BR] habitats dominated by encrusting coralline algae (CCA). Our objective was to assess the ecological functioning of these habitats in the shallow waters of a pristine islet. Our results may then be used to assess any alterations in organic carbon cycles in the coastal benthic communities of these two habitats, in relation to climate change.
section_name: Introduction

section_name: Materials and Methods

section_content: The Dokdo Islets, parts of a volcanic island complex in the East Sea between Japan and Korea, are located near the northeastern boundary of Ulleung Basin (Figure 1 ). They are about 87. 4 km from Ulleung Island and consist of two main islets with submerged volcanic edifices. They rise about 2100 m from the base of Ulleung Basin and have a total area of ~0. 186 km 2 [26]. Unique oceanographic features have been observed around the island. A polar front varies between south and north movement according to the strength of warm and cold waters [27, 28]. Occasionally, eddies formed near the southern Korea peninsular encroach from the east and/or south [29]. These oceanographic dynamics, in the southern part of the basin, may contribute to the high primary production in the region [30]. The macroalgal diversity of Dokdo Islets is greater than that of Ulleung Island and comparable to that of well-developed tidal flats in western Korean waters, reinforcing its status as a biodiversity hotspot [31]. Various seaweed species have been observed on the rocky bottoms of the Dokdo Islets, but these have declined recently, associated with CCA whitening [5, 32, 33]. [30]. The macroalgal diversity of Dokdo Islets is greater than that of Ulleung Isl comparable to that of well-developed tidal flats in western Korean waters, reinfo status as a biodiversity hotspot [31]. Various seaweed species have been observe rocky bottoms of the Dokdo Islets, but these have declined recently, associated w whitening [5, 32, 33].
section_name: Study Sites

section_content: Samplings and AEC experiments were conducted on 10-11 October 2019. Th urements were programmed to start at 07:30 AM (local time) (Table S1 ). The stee depth gradient around the Dokdo Islets creates a narrow zonation of seaweeds an the sites for AEC measurement that need careful selection to represent the benthic properly. With diver assistance, the two measurement sites were chosen using the ing criteria: (1) flat area for firm AEC deployment, (2) dense seaweeds for KB rounded rocks covered with CCA for BR, and (3) representative of the benthic ha
section_name: Site Selection

section_content: The coastal version of the AEC system was developed by the Korea Institute o Science and Technology, KIOST ECI, to estimate O2 exchange rates in benthic bo layers [34]. Briefly, the design was based on the work of Berg et al. [15] for ready tion to coastal-water studies. It consists of four compartments: (1) a frame, (2) dimensional (3D) acoustic Doppler velocimeter (ADV) (Nortek, Vector), (3) a sponse optode sensor (OXR50-UHS, Pyro Science) with an amplifier (Firesting oxygen meter, Pyro Science), and (4) an auxiliary O2 optode sensor for in-situ cal (Anderra, 4330 and 3835) (Figure 2 ). Key modifications of ECI, compared with p AEC systems [15], include a foldable frame and a fast-response O2-sensing syst The unique triangular folding frame enhances portability in the field, thereby a easy operation from a small boat in shallow water. Furthermore, an additional A tem can be installed independently, in parallel, at the opposite side of the AEC in frame [35]. To determine the O2 exchange flux in rocky habitats, high-frequency (<1 Hz) synchronized with 3D velocity measurements are required [15]. The m glass-needle oxygen sensor and amplifier have been widely used for AEC systems However, the durability of the glass sensor limits long-term measurement beca
section_name: Aquatic Eddy Covariance System (AEC)

section_content: Samplings and AEC experiments were conducted on 10-11 October 2019. The measurements were programmed to start at 07:30 AM (local time) (Table S1 ). The steep water depth gradient around the Dokdo Islets creates a narrow zonation of seaweeds and, thus, the sites for AEC measurement that need careful selection to represent the benthic habitats properly. With diver assistance, the two measurement sites were chosen using the following criteria: (1) flat area for firm AEC deployment, (2) dense seaweeds for KB or small rounded rocks covered with CCA for BR, and (3) representative of the benthic habitats.
section_name: Site Selection

section_content: The coastal version of the AEC system was developed by the Korea Institute of Ocean Science and Technology, KIOST ECI, to estimate O 2 exchange rates in benthic boundary layers [34]. Briefly, the design was based on the work of Berg et al. [15] for ready application to coastal-water studies. It consists of four compartments: (1) a frame, (2) a three-dimensional (3D) acoustic Doppler velocimeter (ADV) (Nortek, Vector), (3) a fastresponse optode sensor (OXR50-UHS, Pyro Science) with an amplifier (Firesting O2-Mini oxygen meter, Pyro Science), and (4) an auxiliary O 2 optode sensor for in-situ calibration (Anderra, 4330 and 3835) (Figure 2 ). Key modifications of ECI, compared with previous AEC systems [15], include a foldable frame and a fast-response O 2 -sensing system [35]. The unique triangular folding frame enhances portability in the field, thereby allowing easy operation from a small boat in shallow water. Furthermore, an additional AEC system can be installed independently, in parallel, at the opposite side of the AEC in a single frame [35]. To determine the O 2 exchange flux in rocky habitats, high-frequency O 2 data (<1 Hz) synchronized with 3D velocity measurements are required [15]. The modified glass-needle oxygen sensor and amplifier have been widely used for AEC systems [15, 35]. However, the durability of the glass sensor limits long-term measurement because the sensor commonly breaks and cannot record data for an entire day [35]. To overcome this limitation, the glass oxygen microsensor was replaced with a fast-response (<0. 3 s) optode sensor system [35, 36]. However, the lifetime of the optode and its slow amplifier response time, relative to the glass microsensor system, still limit the AEC system [36]. 

10, x FOR PEER REVIEW 4 of 15 limitation, the glass oxygen microsensor was replaced with a fast-response (<0. 3 s) optode sensor system [35, 36]. However, the lifetime of the optode and its slow amplifier response time, relative to the glass microsensor system, still limit the AEC system [36].
section_name: Aquatic Eddy Covariance System (AEC)

section_content: Each AEC component was assembled as follows: (1) the ADV sensor was aligned vertically, and was either upward-looking in the vegetated habitat (KB) or downwardlooking in the rock bottom (BR) areas; (2) the oxygen sensor tilted at about 45° and its tip located at the edge of the sampling volume of the ADV [36] ; (3) the measurement heights (h) were fixed ~30 cm above the kelp canopy in KB and above the rock surface in BR; (4) the auxiliary O2 optode was independently mounted on the frame at the same measuring height; and (5) the AEC was carefully hand deployed by rope from a small boat with the help of a scuba diver. Using an underwater compass and level, the AEC frame was oriented to the north, without disturbance of the local flow direction, and parallel to the bottom surface. Before AEC device deployment, the ADV was preprogramed, on the boat deck, with a sampling rate of 32 Hz in continuous measuring mode. The photosynthetically available radiation (PAR) at the surface of the habitat was measured independently, at 1 s intervals, using a PAR sensor (JFE Advantech Co, DEFI-L., Nishinomiya, Japan).
section_name: AEC Measurements

section_content: Benthic O2 flux was calculated using custom-developed software written using MATLAB (MathWorks, v. 2018, Natick, MA, USA) and Python (v. 3. 8, Python Software Foundation, Beaverton, OR, USA) following established protocols [36] [37] [38]. Before the derivation of O2 flux, the data were verified and corrected in three sequential steps. In brief, the high-temporal-resolution ADV data were carefully assessed using beam correlations and signal-to-noise ratios (SNRs). Data with poor beam correlation (<70%) and a low SNR (<5) were excluded. The remaining data were then averaged from 32 Hz to 4 Hz to reduce noise and facilitate data handling [36]. Outliers in the downsized 3D velocity and O2 time series data were removed using an acceleration threshold method [39], and these were replaced by linear interpolation between neighboring values. Finally, the filtered data were rotated using a double-rotation method and sliced into 15 min intervals [38]. 

The mean turbulent O2 flux (hereafter O2 flux) was calculated over the 15 min (0. 25 h) periods as:
section_name: Data Analysis and Calculation

section_content: Each AEC component was assembled as follows: (1) the ADV sensor was aligned vertically, and was either upward-looking in the vegetated habitat (KB) or downwardlooking in the rock bottom (BR) areas; (2) the oxygen sensor tilted at about 45 • and its tip located at the edge of the sampling volume of the ADV [36] ; (3) the measurement heights (h) were fixed ~30 cm above the kelp canopy in KB and above the rock surface in BR; (4) the auxiliary O 2 optode was independently mounted on the frame at the same measuring height; and (5) the AEC was carefully hand deployed by rope from a small boat with the help of a scuba diver. Using an underwater compass and level, the AEC frame was oriented to the north, without disturbance of the local flow direction, and parallel to the bottom surface. Before AEC device deployment, the ADV was preprogramed, on the boat deck, with a sampling rate of 32 Hz in continuous measuring mode. The photosynthetically available radiation (PAR) at the surface of the habitat was measured independently, at 1 s intervals, using a PAR sensor (JFE Advantech Co., DEFI-L., Nishinomiya, Japan).
section_name: AEC Measurements

section_content: Benthic O 2 flux was calculated using custom-developed software written using MAT-LAB (MathWorks, v. 2018, Natick, MA, USA) and Python (v. 3. 8, Python Software Foundation, Beaverton, OR, USA) following established protocols [36] [37] [38]. Before the derivation of O 2 flux, the data were verified and corrected in three sequential steps. In brief, the high-temporal-resolution ADV data were carefully assessed using beam correlations and signal-to-noise ratios (SNRs). Data with poor beam correlation (<70%) and a low SNR (<5) were excluded. The remaining data were then averaged from 32 Hz to 4 Hz to reduce noise and facilitate data handling [36]. Outliers in the downsized 3D velocity and O 2 time series data were removed using an acceleration threshold method [39], and these were replaced by linear interpolation between neighboring values. Finally, the filtered data were rotated using a double-rotation method and sliced into 15 min intervals [38]. 

The mean turbulent O 2 flux (hereafter O 2 flux) was calculated over the 15 min (0. 25 h) periods as:

where the overbars indicate time averaging and w and c represent the instantaneous fluctuations away from the mean of the vertical velocity (cm s -1 ) and O 2 concentration (mmol m -3 ). The fluctuation components, w and c, were extracted using Reynolds decomposition, where c = cc and w = ww with w and c being, respectively, the measured vertical velocity and oxygen in situ, and w and c being the time averaged values [15]. The decomposition was performed by linear detrending using MATLAB's detrend function. The window size of 15 min has been defined in previous studies as the optimal time period for maintaining a constant flux signal [17]. Then, hourly mean O 2 flux was calculated by averaging four 15 min O 2 fluxes [16]. 

Categories of light and dark were grouped using averaged 15 min PAR values while assuming that the critical value of PAR across a day was >0. 2 µmol photon m -2 d -1. The benthic respiration, GPP, and NEM were calculated for each 15 min window as follows: R = 1 n

where n is the number of 15 min fluxes, ∑ flux dark is the flux during the nighttime (mmol O 2 m -2 d -1 ), ∑ flux light is the flux during the daytime (mmol O 2 m -2 d -1 ), h light is the number of hours of light, and h dark is the number of hours of darkness [16]. 

Assuming that the O 2 flux in the benthic habitats was the product of photosynthesis (P) and respiration (R) by benthic organisms, the light dependencies of the O 2 fluxes (P-I curves) were estimated at both sites using the relationship between the hourly mean daytime O 2 flux and PAR [40]. Specifically:

where P max is the maximum photosynthesis rate (mmol O 2 m -2 d -1 ), I is the measured PAR at the canopy height in KB or at the bottom surface in BR (µmol photon m -2 s -1 ), I k is the light saturation constant (µmol photon m -2 s -1 ), and R d is the dark respiration (mmol O 2 m -2 d -1 ). The hourly mean O 2 fluxes were fitted to the 1 h mean PAR. The best fitting parameters, P max, I k, and R d, were estimated using the "curve_fit" function of the scipy module in Python. Additional parameters, α (photosynthetic efficiency) and I c (compensation irradiance), were calculated as P max /I k and R d /α, respectively.
section_name: Data Analysis and Calculation

section_content: To determine the species richness and biomass of macrophytes in KB, a 50 × 50 cm square frame was installed randomly at two spots for duplication. Samples, including epiphytes, on the rock surface in the square frame were carefully collected and fixed using natural formalin. The samples were transported to the laboratory and the macrophytes were identified to species using a microscope. Their wet biomass was measured with a balance, and the number of individuals was counted using image analysis.
section_name: Sampling and Species Richness

section_name: Results

section_content: In KB, the dominant kelp species was Eisenia bicyclis, which had wet biomass of 2710 g m -2. This was lower than previously measured at Dokdo (2000: 7600-29,000 g m -2, 2014: 895-5377 g m -2 ) [5]. Eisenia bicyclis individuals ranged from 357 to 1656 g, with an average of 678 g wet weight. The total length of the thallus ranged from 43. 0 to 85. 0 cm, with an average of 60. 2 cm. The typical CCA species observed at the bottom in KB were Synarthrophyton chejuensis and Lithophyllum okamurae [6, 33]. There were no macroalgae in BR, but CCA were abundant on the rock surface. The mollusks (snails and mussels) and starfish were sparsely distributed. 

Mean KB water depth was 5. 1 m, which was shallower than BR (8. 4 m) (Table S1 ). The bottom water temperature varied from 18. 83 to 19. 63 • C (mean ± standard deviation [SD], 19. 20 ± 0. 20 • C) and salinity from 33. 82 to 34. 26 (mean ± SD, 34. 02 ± 0. 07) (Table S1 ). The daytime PAR was up to 493 µmol photon m -2 s -1 (mean: 74 µmol photon m -2 s -1 ) for KB and up to 399 µmol photon m -2 s -1 (average: 67 µmol photon m -2 s -1 ) for BR, which might depend on water depth.
section_name: Site Characteristics

section_content: Using the continuous measuring mode of the ADV, 3D velocity and O 2 were measured for ~25 h (Table S1 ). To avoid the inclusion of abnormal data during deployment and recovery of the AEC, both ends of the data were excluded from the analysis; thus, the measurements used for further calculation were less than 24 h. 

As expected, the 3D velocity was significantly higher in KB than in BR as a consequence of wave orbital motion (Figure 3 ) [35, 41, 42]. The mean flow velocity (U = x 2 + y 2 + z 2 ) in KB was estimated as 8. 8 cm s -1, which was about twice that in BR (3. 6 cm s -1 ) (Table 1 ). S1 ). 

The daytime PAR was up to 493 μmol photon m -2 s -1 (mean: 74 μmol photon m -2 s -1 ) for KB and up to 399 μmol photon m -2 s -1 (average: 67 μmol photon m -2 s -1 ) for BR, which might depend on water depth.
section_name: Time Series Data

section_content: Using the continuous measuring mode of the ADV, 3D velocity and O2 were measured for ~25 h (Table S1 ). To avoid the inclusion of abnormal data during deployment and recovery of the AEC, both ends of the data were excluded from the analysis; thus, the measurements used for further calculation were less than 24 h. 

As expected, the 3D velocity was significantly higher in KB than in BR as a consequence of wave orbital motion (Figure 3 ) [35, 41, 42]. The mean flow velocity ( = + + ) in KB was estimated as 8. 8 cm s -1, which was about twice that in BR (3. 6 cm s -1 ) (Table 1 ). Oxygen time series did not show any temporal variation at either site, but the mean O2 concentration was significantly higher (t-test, p < 0. 05) in KB (277 μmol L -1 ) than in BR (227 μmol L -1 ) (Figure 4, Table S1 ), which may be an effect of O2 production by benthic flora (mainly E. cava and E. bicyclis) photosynthesis. Abnormal spikes in the O2 time series data were observed around 12-13 h in KB, which may be interpreted as interference by debris particles. Similar random spikes in BR were excluded from further analysis (Figure 4 ). Oxygen time series did not show any temporal variation at either site, but the mean O 2 concentration was significantly higher (t-test, p < 0. 05) in KB (277 µmol L -1 ) than in BR (227 µmol L -1 ) (Figure 4, Table S1 ), which may be an effect of O 2 production by benthic flora (mainly E. cava and E. bicyclis) photosynthesis. Abnormal spikes in the O 2 time series data were observed around 12-13 h in KB, which may be interpreted as interference by debris particles. Similar random spikes in BR were excluded from further analysis (Figure 4 ). The diurnal cumulative O2 flux showed overall opposite trends during daytime and nighttime (Figure 5 ). The cumulative O2 flux increased during the daytime because of consistent O2 production by benthic flora photosynthesis, but decreased at nighttime because of respiration of benthic organisms. Cumulative fluxes varied with PAR intensities and the magnitude of variation was significantly higher in KB than in BR.
section_name: Time Series Data

section_content: Hourly mean O2 fluxes were clearly related to diurnal PAR changes (Figure 6 ). The hourly mean O2 fluxes ranged from -245 to 475 mmol O2 m -2 d -1 (mean ± standard error (SE): 2 ± 42 mmol O2 m -2 d -1 ) in KB and from -95 to 172 mmol O2 m -2 d -1 (mean ± SE: -14 ± 20 mmol O2 m -2 d -1 ) in BR. Differences in fluxes between sites (KB vs. BR) and times (daytime vs. nighttime) were statistically significant (t-test, p < 0. 01), suggesting that the pho- The diurnal cumulative O2 flux showed overall opposite trends during daytime and nighttime (Figure 5 ). The cumulative O2 flux increased during the daytime because of consistent O2 production by benthic flora photosynthesis, but decreased at nighttime because of respiration of benthic organisms. Cumulative fluxes varied with PAR intensities and the magnitude of variation was significantly higher in KB than in BR.
section_name: Hourly Averaged O2 Flux with PAR

section_content: Hourly mean O2 fluxes were clearly related to diurnal PAR changes (Figure 6 ). The relationship between O2 flux and PAR is shown in Figure 7. The parameters of the in-situ P-I curve are listed in Table 1. The modeled maximum net O2 production (Pmax) and dark respiration (Rd) were 50% higher in KB than in BR, which implies that both photosynthesis and respiration are higher in the kelp bed than in the CCA mat. Furthermore, the light saturation (Ik) was significantly lower in KB than in BR, which suggests that kelp responds to irradiance more efficiently. The relationship between O 2 flux and PAR is shown in Figure 7. The parameters of the in-situ P-I curve are listed in Table 1. The modeled maximum net O 2 production (P max ) and dark respiration (R d ) were 50% higher in KB than in BR, which implies that both photosynthesis and respiration are higher in the kelp bed than in the CCA mat. Furthermore, the light saturation (I k ) was significantly lower in KB than in BR, which suggests that kelp responds to irradiance more efficiently. The relationship between O2 flux and PAR is shown in Figure 7. The param the in-situ P-I curve are listed in Table 1. The modeled maximum net O2 productio and dark respiration (Rd) were 50% higher in KB than in BR, which implies that bo tosynthesis and respiration are higher in the kelp bed than in the CCA mat. Furthe the light saturation (Ik) was significantly lower in KB than in BR, which suggests th responds to irradiance more efficiently.
section_name: Hourly Averaged O2 Flux with PAR

section_name: Pmax

section_content: The relationships between O 2 fluxes, which were binned in 1 cm s -1 intervals for clarity, and mean velocities (U) per bin in KB and BR, are shown in Figure 8. The variabilities in the O 2 fluxes and U values were larger in KB than in BR. Overall, the O 2 fluxes were positive during daytime and negative at nighttime, which may be the result of the benthic biological activity (photosynthesis and respiration). Apart from in the daytime in KB, the relationships between O 2 flux and U were linear, increasing with U in daytime, but decreasing with U at night. Mar. Sci. Eng. 2022, 10, x FOR PEER REVIEW biological activity (photosynthesis and respiration). Apart from in the dayti relationships between O2 flux and were linear, increasing with in day creasing with at night.
section_name: O 2 Flux and Mean Current

section_content: Considerable differences in GPP and R were observed between the KB tats (Figure 9, Table 1 ). The GPP in KB was 163 mmol O2 m -2 d -1, which wa times higher than that in BR (51 mmol O2 m -2 d -1 ). Similarly, R was higher in K O2 m -2 d -1 ) than in BR (63 mmol O2 m -2 d -1 ). The NEMs for KB and BR were 8 mmol O2 m -2 d -1 and -12 mmol O2 m -2 d -1, suggesting that KB is net autotro BR is net heterotrophic.
section_name: Net Ecosystem Metabolism (NEM)

section_content: Considerable differences in GPP and R were observed between the KB and BR habitats (Figure 9, Table 1 ). The GPP in KB was 163 mmol O 2 m -2 d -1, which was about three-times higher than that in BR (51 mmol O 2 m -2 d -1 ). Similarly, R was higher in KB (154 mmol O 2 m -2 d -1 ) than in BR (63 mmol O 2 m -2 d -1 ). The NEMs for KB and BR were, respectively, 8 mmol O 2 m -2 d -1 and -12 mmol O 2 m -2 d -1, suggesting that KB is net autotrophic, whereas BR is net heterotrophic. 

times higher than that in BR (51 mmol O2 m -2 d -1 ). Similarly, R was high O2 m -2 d -1 ) than in BR (63 mmol O2 m -2 d -1 ). The NEMs for KB and BR 8 mmol O2 m -2 d -1 and -12 mmol O2 m -2 d -1, suggesting that KB is net au BR is net heterotrophic.
section_name: Net Ecosystem Metabolism (NEM)

section_name: Discussion

section_content: Our data clearly indicated that O 2 fluxes varied with flow velocity at both sites (Figure 8 ). The thickness of the diffusive boundary layer at the solid-water interface depends on the friction force: an increasing flow rate can narrow the diffusion boundary layer, which can stimulate both production and respiration in benthic communities [10, 16, [19] [20] [21]. Indeed, studies on the metabolism of seagrass (eelgrass) meadows show that local flow velocity can influence O 2 flux as much as six-fold [16]. Such results are consistent across benthic habitats [17, 22, 23]. At night, the influence of flow on O 2 flux was about five-times higher in KB than in BR, which suggests that wave forcing, caused by wave orbitals, facilitates kelp respiration [20, 43, 44].
section_name: O 2 Flux Variability

section_content: At Dokdo, the two studied habitats had quite different hyperbolic P-I relationships (Figure 7 ). P max and R d were significantly higher in KB than in BR (Table 1 ), whereas the I k and I c in KB were each about one-third of their values in BR. This implies that photosynthesis is more responsive to light in kelp than in CCA. 

Murakami et al. [45] evaluated the seasonal photosynthetic levels of E. cava and E. bicyclis through laboratory incubation experiments. They found little seasonal change in photosynthetic parameters (P max, R d, I k, and I c ). P max ranged from 257 to 321 mmol O 2 m -2 d -1 and R d from 62 to 73 mmol O 2 m -2 d -1, values that are similar to our findings (P max : 257 mmol O 2 m -2 d -1, R d : 65 mmol O 2 m -2 d -1 ) (Table 1 ). Their net photosynthesis (O 2 flux) was saturated at light levels above 100 µmol photon m -2 s -1, which closely matched our own in-situ measurements. However, their light saturation level parameters (I k and I c ) were 30%-40% below our results (Table 2 ). In a previous study, the I k and I c of the leaf or thallus tissues of selected plants were significantly lower than overall community values, which was attributed to internal self-shading in the community [46]. To date, incubation experiments for the estimation of the photosynthesis of aquatic macrophytes have used an enclosed light (-dark) chamber under an irradiance gradient, either in the laboratory or in situ [11, 12, 45, 46]. As the community-level characteristics of photosynthesis are tightly coupled with kelp density, light absorption, and light distribution within the bed, incubation results using phytoelements may need to be more carefully considered to understand benthic photosynthesis at the community scale [46]. On the other hand, AEC results can noninvasively assess integrated O 2 dynamics in whole benthic communities on a large footprint (10 to 50 m 2 ) [25]. AEC application may, therefore, provide valuable information about the factors influencing photosynthesis in benthic habitats over multiple time scales (up to yearly) [19].
section_name: P-I Relationship

section_content: The NEMs of macroalgal communities in the Northern Hemisphere are listed in Table 2 [ Our results are comparable with other macroalgal communities, even though taken in a low-growth period. 

The NEM of a specific target community is an important indicator of local ecological function [51]. The community-wide difference between GPP and R (NEM = GPP -R) can be a proxy for carbon flow in the food web [51]. For example, a positive NEM indicates a net autotrophic condition, and reflects the degree of organic matter accumulation within the community and, hence, the excess organic carbon available to nearby ecosystems as food. By contrast, a negative NEM signals net heterotrophy and indicates the degree to which a community requires organic carbon, from an exterior ecosystem, as food-which also implies that it acts as a carbon sink. Our results clearly indicate that organic carbon is produced in kelp beds, and some of this is exported to areas with little or no vegetation. This suggests a "production-consumption" relationship between benthic communities in a localized area. Our methods and results do not, however, permit calculation of carbon flow between KB and BR.
section_name: NEM of Dokdo Islets

section_content: Kelp beds in shallow coastal waters are among the most productive ecosystems in the world [52], and serve as nurseries for juveniles, an energy source for the food web, and a carbon source for adjacent habitats [53] [54] [55] [56] [57] [58]. However, kelp beds are being progressively destroyed by disease, herbivores (mainly sea urchins), and physiological stresses from various environmental factors (e. g., suboptimal light, temperature, nutrients, substrate), including climate change [1, 52]. In Japan, the periodic intrusion of the warm Tsushima Current along the coast produces mass mortality of Ecklonia and Eisenia species in rocky habitats, which are replaced by tropical species (Sargassum ilicifolium) on the small spatial scale for a short time [6]. Similar to the situation in Japan, kelp beds are rapidly disappearing from the coastlines of the Korean peninsula and are replaced by the whitening of CCA [7, 59, 60]. Considering that the increase in the annual mean surface water temperature around Korea is about three-times higher than that of the global trend over the last five decades, the physiological interaction in the benthic community by thermal effect may be closely related to the contraction of kelp and expansion of CCA [61]. 

These changes bring risks to entire benthic habitats and ecosystem structures. Although Dokdo is less impacted by anthropogenic activity because of its geological setting, changes in oceanographic characteristics under climate change can threaten its benthic ecosystems [1, 7, 8]. Thus, the assessment of community metabolisms in kelp beds and whitening areas is important to gain a better understanding of their ecological functioning and, hence, to aid in the conservation of pristine conditions at Dokdo and, more broadly, the understanding of carbon cycles around open ocean islets. 

The Dokdo Islets remain relatively pristine because of long-term limited human access and are, therefore, an excellent location to monitor environmental changes. Recently, the coastal ecosystems in the north Pacific have been altered by climate change [1, 56], but quantitative assessment of the changes remains limited in benthic communities. Our results suggest that the trophic status of rocky habitats can be determined from their vegetation, which in the case of kelp beds, can act as both a carbon source and sink over small horizontal scales. The present study was implemented when the growth of macroalgae was low, and, thus, further long-term research is needed to assess fully the ecological functioning of rocky habitats and the biogeochemical cycles of organic carbon in coastal waters during different seasons and across the year. In addition, work such as that described in this paper needs to be undertaken at replicate sites and with replication within habitats, for statistically robust inferences to be drawn.
section_name: Implications and Further Suggestions

section_content: In the present study, for the first time, in-situ AEC measurements were used to estimate the NEM, GPP, and R in the typical benthic habitats (kelp bed and bare rocky bottom) of the Dokdo Islets. In-situ P-I relationships were also assessed to gain a better understanding of the major factors controlling organic carbon cycles in the euphotic waters of the islets. Our results indicate that O 2 flux variations, over multiple time scales (from 15 min to 1 day), are driven by local flow rates and irradiance. The O 2 fluxes have a strong positive linear relationship with local flow. The P-I curves indicate that the rates of photosynthetic activities are significantly higher in the kelp beds (KB) than in the relatively unvegetated bottom habitats (BR) dominated by CCA. The photosynthesis light response of kelp was more efficient than that of CCA. The NEM of the benthic communities at Dokdo may be determined primarily by the nature of the vegetation, as CCA-dominated habitats have consistently reduced metabolism parameters relative to kelp beds. The trophic states were assessed as net autotrophic on the kelp bed, but net heterotrophic on the bare rocky bottom. Application of these in-situ AEC methods can be expected to provide robust benthic metabolism rates for kelp beds, as well as for CCA mats in hard-bottom habitats.
section_name: Conclusions

section_name: Acknowledgments:

section_content: We wish to thank Min for assistance during AEC deployment and seaweed sampling. The support from Park and Lee was helpful and greatly appreciated. In addition, we also thank the three anonymous reviewers for their helpful comments that improved the early version of the MS.
section_name: Acknowledgments:

section_content: The authors declare no conflict of interest.
section_name: Conflicts of Interest:

section_content: Funding: This work was supported by the Korea Institute of Ocean Science and Technology ( PEA0012 ), the Ministry of Oceans and Fisheries ("a sustainable research and development of Dokdo", PG52911 ), and the National Institute of Fisheries Science ( R2022062 ). 

Institutional Review Board Statement: Not applicable.

section_content: The following supporting information can be downloaded at: https: //www. mdpi. com/article/10. 3390/jmse10070887/s1, Table S1 : AEC measurement times and

section_name: Data Availability Statement:

section_content: The data used to support the findings of this study are available from the corresponding author upon request.
section_name: Data Availability Statement:

section_content: Informed Consent Statement: Not applicable.