Integrating Aquatic Metabolism and Net Ecosystem CO2 Balance in Short- and Long-Hydroperiod Subtropical Freshwater Wetlands

How aquatic primary productivity influences the carbon (C) sequestering capacity of wetlands is uncertain. We evaluated the magnitude and variability in aquatic C dynamics and compared them to net ecosystem CO2 exchange (NEE) and ecosystem respiration (Reco) rates within calcareous freshwater wetlands in Everglades National Park. We continuously recorded 30-min measurements of dissolved oxygen (DO), water level, water temperature (Twater), and photosynthetically active radiation (PAR). These measurements were coupled with ecosystem CO2 fluxes over 5 years (2012–2016) in a long-hydroperiod peat-rich, freshwater marsh and a short-hydroperiod, freshwater marl prairie. Daily net aquatic primary productivity (NAPP) rates indicated both wetlands were generally net heterotrophic. Gross aquatic primary productivity (GAPP) ranged from 0 to − 6.3 g C m−2 day−1 and aquatic respiration (RAq) from 0 to 6.13 g C m−2 day−1. Nonlinear interactions between water level, Twater, and GAPP and RAq resulted in high variability in NAPP that contributed to NEE. Net aquatic primary productivity accounted for 4–5% of the deviance explained in NEE rates. With respect to the flux magnitude, daily NAPP was a greater proportion of daily NEE at the long-hydroperiod site (mean = 95%) compared to the short-hydroperiod site (mean = 64%). Although we have confirmed the significant contribution of NAPP to NEE in both long- and short-hydroperiod freshwater wetlands, the decoupling of the aquatic and ecosystem fluxes could largely depend on emergent vegetation, the carbonate cycle, and the lateral C flux.

ABSTRACT How aquatic primary productivity influences the carbon (C) sequestering capacity of wetlands is uncertain. We evaluated the magnitude and variability in aquatic C dynamics and compared them to net ecosystem CO 2 exchange (NEE) and ecosystem respiration (R eco ) rates within calcareous freshwater wetlands in Everglades National Park. We continuously recorded 30-min measurements of dissolved oxygen (DO), water level, water temperature (T water ), and photosynthetically active radiation (PAR). These measurements were coupled with ecosystem CO 2 fluxes over 5 years (2012)(2013)(2014)(2015)(2016) in a long-hydroperiod peat-rich, freshwater marsh and a short-hydroperiod, freshwater marl prairie. Daily net aquatic primary productivity (NAPP) rates indicated both wetlands were generally net heterotrophic. Gross aquatic primary productivity (GAPP) ranged from 0 to -6.3 g C m -2 day -1 and aquatic respiration (R Aq ) from 0 to 6.13 g C m -2 day -1 . Nonlinear interactions between water level, T water , and GAPP and R Aq resulted in high variability in NAPP that contributed to NEE. Net aquatic primary productivity accounted for 4-5% of the deviance explained in NEE rates. With respect to the flux magnitude, daily NAPP was a greater proportion of daily NEE at the long-hydroperiod site (mean = 95%) compared to the short-hydroperiod site (mean = 64%). Although we have confirmed the significant contribution of NAPP to NEE in both long-and shorthydroperiod freshwater wetlands, the decoupling of the aquatic and ecosystem fluxes could largely depend on emergent vegetation, the carbonate cycle, and the lateral C flux.

INTRODUCTION
Wetland ecosystems worldwide are recognized for their disproportionately large and complex organic carbon (C) storage potential (Mitsch and Gosselink 2007;Nahlik and Fennessy 2016;Lu and others 2017). Although they are believed to be highly productive ecosystems with low rates of decomposition due to anaerobic conditions (Chmura and others 2003), wetlands vary in their capacity to sequester C based on hydrology (Kayranli and others 2010;Bernal and Mitsch 2012;Malone and others 2014;Lu and others 2017). Although eddy covariance and static chamber studies have been conducted to understand ecosystem variability and the relative contributions of different land cover types to C sequestration (Morin and others 2017; Rey-Sanchez and others 2018), it is not clear how C pools and fluxes from emergent and submerged parts of wetland ecosystems interact to influence net C storage.
An indicator of the current potential for C sequestration, net ecosystem exchange of CO 2 (NEE) is the balance between gross ecosystem exchange (GEE) of C and ecosystem respiration (R eco ). In a global synthesis of wetlands, Lu and others (2017) showed that the leaf area index of emergent vegetation was correlated with GEE and R eco , but not with NEE. This suggests that although the GEE and R eco of wetlands is heavily influenced by emergent vegetation (Lu and others 2017), the hydrological control on the balance between the two flux components may be the dominant driver of the C sequestering capacity of wetland ecosystems (Bernal and Mitsch 2012). In addition to creating anaerobic conditions that slow decomposition and C loss, inundation also supports C uptake through aquatic primary productivity. Net aquatic primary productivity (NAPP) specifically refers to the balance between C uptake by periphyton, submerged aquatic plants, algae, and cyanobacteria, and the respiration of aquatic vegetation (R Aq ) and decomposition (Davis and Ogden 1994;McCormick and Laing 2003). Depending on the dominance of gross aquatic primary productivity (GAPP) or R Aq , the NAPP can enhance the sink or source capacity of ecosystems. However, there have been few studies that have shown how aquatic primary productivity influences net fluxes of CO 2 from wetlands at the ecosystem scale (but see Hagerthey and others 2010).
Additionally, the impact of calcium carbonate (CaCO 3 ) production and dissolution on the C cycle is not well-understood, further complicating C dynamics in calcareous wetland ecosystems. Calcium carbonate production in freshwater ecosystems is partially attributed to calcifying algae (Wefer 1980;Rupp and Adams 1981;Heath and others 1995;Davis and others 2005). Although calcification sequesters C, it is also a source of CO 2 to the surrounding water (Frankignoulle and others 1994). A decline in total alkalinity facilitates the return of CO 2 to the water column (Frankignoulle and others 1994). There is still much uncertainty about the fate of CO 2 released after calcification, as well as how much of this CO 2 is used in GAPP, or is respired (Macreadie and others 2017).
Hydroperiod is a critical driver of C storage in wetlands, but the relative contributions of water column metabolism to net ecosystem C storage are not well-constrained. The complexity of wetland CO 2 dynamics requires multiple ecosystem-level measurements from both within and above the water column to adequately understand the C sequestering capacity of these ecosystems. As a result, few studies (Hagerthey and others 2010) have evaluated whole-wetland ecosystem CO 2 fluxes to understand the role of aquatic CO 2 dynamics on the C sequestering capacity of wetland ecosystems. In the subtropical Florida Everglades, NEE varies across wetland types and in response to shifts in the physiological activity of emergent vegetation with seasonal changes in water level others 2010, 2012;Jimenez and others 2012;Malone and others 2014). These calcareous freshwater wetlands fluctuate from being a small source to a small sink of CO 2 annually, which is strongly linked to seasonal subtropical hydrology (Schedlbauer and others 2010, 2012; Jimenez and others 2012; Malone and others 2014; Zhao and others 2019). Seasonal changes in photosynthetic rates of emergent plant species partially explain patterns in NEE (Zhao and others 2018), yet metabolism within the water column may also have an effect on how much and how long C is sequestered. Although NAPP in the Everglades has been studied in response to variation in phosphorus concentrations (Davis and Ogden 1994;McCormick and Laing 2003;Iwaniec and others 2006) and habitat type (Ewe and others 2006;Hagerthey and others 2010), the magnitude and drivers of temporal variation in aquatic C dynamics in freshwater wetlands remains uncertain. Further, the Florida Everglades have been profoundly transformed through water diversion and land conversion over the last century (Light and Dineen 1994;Troxler and others 2013), adding to the complexity of predicting future C storage in wetlands.
In response to wetland degradation, current water management practices are being modified under the Comprehensive Everglades Restoration Plan (CERP), with the goal to re-establish hydroperiods closer to natural seasonal regimes and reduce chronically low water levels (Perry 2004). These hydrologic changes introduce more complexity and uncertainty into our understanding of the system, making it essential that we determine how activity within the water column will change with these alterations and their effects on biogeochemical cycles (Davis and Ogden 1994). The ongoing landscape-scale restoration of the Everglades (Sklar and others 2005) is expected to increase freshwater flow, but the degree to which this will change the C sequestration capacity is unknown (Troxler and others 2013).
Here, we evaluated the magnitude and variability in aquatic C dynamics and compared this to NEE and R eco in Everglades calcareous freshwater wetlands. In these wetlands where hydroperiods are long enough to limit the breakdown of submerged organic material, C is mainly sequestered into peat. (Davis and Ogden 1994;Gleason and Stone 1994;Davis and others 2005). In short-hydroperiod marl prairie ecosystems, flocculent organic matter (floc), which is formed from periphyton and submerged aquatic vegetation, is highly labile (Pisani and others 2013), resulting in rapid organic matter breakdown during the dry periods and the development of marl soils (Chen and others 2013). We hypothesized that: 1) within Everglades calcareous freshwater wetlands, aquatic primary productivity would increase with increasing water level and air temperature, and 2) there would be weaker relationships between NEE and NAPP in the marl prairie compared to the freshwater marsh where more CO 2 is produced via calcification, which can be used for aquatic photosynthesis in the water column or released to the atmosphere.
In this subtropical region, freshwater wetlands have a year-round growing season. Most precipitation falls in the wet season, which extends from May to October. However, seasonal rainfall patterns combined with elevation, location, and water management practices generate variation in water levels and the length of annual inundation periods.

Dissolved Oxygen and Ancillary Micrometeorological Measurements
Each site was instrumented with a D-Opto SDI-12 Optical DO sensor (ENVCO Global, New Zealand) in 2011 to determine dissolved oxygen (DO) concentration (ppm) via its relationship with fluorescence of ruthenium and temperature in the water column. At the long-hydroperiod marsh site (SRS-2), the DO sensor was placed in the slough, mounted about 20 cm from the soil surface on the end of a boardwalk. At the marl Integrating Aquatic Metabolism and Net Ecosystem prairie site (TS/Ph-1), the sensor was placed in a shallow solution pond about 20 cm from the bottom that allowed the sensor to remain submerged even when the water level dropped below the soil surface; however, data were only used when the site as a whole was inundated, allowing continuous mixing with the surrounding water. Although there was no sawgrass and muhly grass within the solution hole, periphyton and aquatic plants were abundant, and the O 2 concentrations in the flowing water partly reflect metabolism of vegetation that is upstream of the sensor. It is important to note that this design ensures changes in O 2 concentrations are reflecting aquatic activity.

Aquatic Ecosystem Metabolism
Prior to calculating NAPP, all 30-min DO data were QA/QC and filtered to remove outliers. The availability of 30-min DO data differed between both sites (we filtered out 29% at SRS-2 and 25% at TS/ Ph-1 of the 5 years) due to equipment outages and weather events which prevented site access. An additional 32% of DO data were removed from the time series at TS/Ph-1 due to low water levels. The continuously measured diel changes in DO concentrations at 30-min intervals were used to estimate NAPP (Odum 1956;Cole and others 2000;Huggins 2005; Hagerthey and others 2010; Argerich and others 2016). Changes in DO are a result of two processes, NAPP and diffusive exchange with the atmosphere or reaeration (D) (Hagerthey and others 2010).
where DDO is the change in measured oxygen concentration from one half hour to the next, d is the water depth, O 2 sat is the theoretical oxygen concentration in the water at equilibrium with the atmosphere estimated as a function of atmospheric pressure and air temperature (Garcia and Gordon 1992). The value K is the coefficient of gas exchange for O 2 at a given temperature, and in aquatic systems, it is generally modeled as a function of wind speed (w) (Cole and others 2000), although studies have shown that commonly used wind speed/gas exchange relationships may overestimate the gas transfer velocity in wetland ecosystems (Ho and others 2018). To estimate GAPP and R A , we fit a model of diel metabolism to the measured DO data using the nls function in R (Hall and others 2015), utilizing the nonlinear model of Van de Bogert and others (2007) with 30min data: where O 2i is the O 2 at time i, O 2i-1 is the O 2 at the previous time, z is the water level (m), and Dt i is the measurement interval of logged O 2 data. In this model, the estimated coefficient R Aq is a negative O 2 flux because O 2 is being consumed, while the estimated coefficient GAPP is a positive O 2 flux. Major assumptions of this approach are that GAPP is a linear function of light (Van de Bogert and others 2007) and that R Aq is constant throughout the day. The model was fit using days where 100% of the data were non-missing (SRS-2, n = 725 days; TS/Ph-1, n = 470 days) to estimate the daily rate of GAPP and R Aq . This resulted in models with the correct likelihood, where the asymptotic likelihood inference was valid. The estimates of GAPP, R Aq , and NAPP were then changed to atmospheric convention (R Aq positive and GAPP negative) to facilitate comparisons with the ecosystem exchange rates. To compare aquatic metabolism to ecosystem metabolism, daily rates of NAPP, R Aq , and GAPP were converted to g C m -2 day -1 assuming a C/O 2 molar ratio of 0.375 and a photosynthetic quotient of 1.2 (McCormick and others 1997; Iwaniec and others 2006; Hagerthey and others 2011).

Ecosystem CO 2 Dynamics
At the ecosystem scale, NEE provides a view of the CO 2 sequestering capacity of wetlands, integrating the exchange of CO 2 with the atmosphere for both organic and inorganic processes. Eddy covariance towers at each site were instrumented with an open-path infrared gas analyzer (IRGA, LI-7500, Li-COR Inc., Lincoln, NE) to measure CO 2 concentration (C; mg mol -1 ) and water vapor molar density (qv; mg mol -1 ), and a paired sonic anemometer (CSAT3, Campbell Scientific Inc., Logan, UT) to measure sonic temperature (T s ; K) and three-dimensional wind speed (u, v and w, respectively, m s -1 ). These sensors were 0.09 m apart and installed at 3.30 and 3.24 m above ground level (a.g.l.) at TS/Ph-1 and SRS, respectively. Data were logged at 10 Hz on a datalogger (CR1000, Campbell Scientific Inc.) and stored on 2 GB CompactFlash cards. The LI-7500's were calibrated monthly using a trace gas standard for CO 2 in air (+ 1.0%), dry N 2 gas and a portable dew point generator (LI-610, LI-COR Inc.). Footprint analyses others 2002, 2004) indicated that 80% of measured fluxes were within 100 m of the tower during convective conditions at both sites. Other meteorological variables were measured at 1-s and collected as half-hourly means, acquired by the same datalogger, and included air temperature (T air ;°C) and relative humidity (RH; %) (HMP45C, Vaisala, Helsinki, Finland) mounted within an aspirated shield (43,502, R.M. Young Co., Traverse City, MI), and barometric pressure (P; atm) (PTB110, Vaisala). The T air /RH sensors were installed at the same height a.g.l. as the IRGA and CSAT.
From 2008 to 2016, 10 Hz raw flux data were processed with EdiRe (v. 1.4.3.1184) (Clement 1999), which included despiking and air density corrections (Webb and others 1980;Aubinet and others 1999). Fluxes (NEE, H, LE) were then corrected for mass transfer resulting from changes in density not accounted for by the IRGA (Webb and others 1980;Massman 2004), and barometric pressure data were used to correct the fluxes to standard atmospheric pressure.
All measurements were filtered when systematic errors in fluxes were indicated, such as (1) evidence of rainfall, condensation, or bird fouling in the sampling path of the IRGA or sonic anemometer, (2) incomplete 30-min datasets during system calibration or maintenance, (3) poor coupling of the canopy with the external atmospheric conditions, as defined by the friction velocity, u*, using a threshold, 0.15 m s -1 (Goulden and others 1996; Clark and others 1999), and (4) excessive variation from the half-hourly mean based on an analysis of standard deviations for u, v, and w wind and CO 2 statistics. Quality assurance of the flux data was also maintained by examining plausibility tests for implausible NEE values (I.e., NEE < -30 or > 30 lmol m -2 s -1 ), as previous research has shown such values are beyond the limits of physiological activity in these wetlands (Jimenez and others 2012; Malone and others Integrating Aquatic Metabolism and Net Ecosystem 2014), as well as stationarity criteria and integral turbulent statistics (Foken and Wichura 1996;Foken and Leclerc 2004). At SRS, 20% of daytime data and 59% of nighttime data were removed. At TS/Ph-1, equipment outages and power failures led to higher amounts of missing or removed data (61% and 75% of the day and nighttime data, respectively).
Following (Jimenez and others 2012;others 2014, 2016)), missing 30-min NEE data were gap-filled using separate functions for day and night. During the day, NEE was gap-filled using a Michaelis-Menten saturation model (NEE day ; Eq. 5) (Michaelis and Menten 1913;Hollinger and others 1999;Falge and others 2001;Ruppert and others 2006), and at night NEE was gap-filled using an Arrhenius equation (NEE night ; Eq. 6): where a is the apparent quantum efficiency, u is PAR, R eco is ecosystem respiration (lmol CO 2 m -2 s -1 ), P max is the maximum ecosystem CO 2 uptake rate (lmol CO 2 m -2 s -1 ), R 0 is the base respiration rate when the air temperature is 0°C and b is an empirical coefficient. In Eq. 5, R eco is an estimated model parameter, whereas R eco measurements are the dependent variable in Eq. 6. Equations were fit monthly using the nls function in R. In cases where monthly data were missing, seasonal (within and between years) or annual models were used (Jimenez and others 2012; Malone and others 2014). Differences in timescales (monthly, seasonally, or annually) for model fit allow model parameters to vary with conditions. Following gap filling, R eco during the day was estimated using Eq. 6, and GEE was calculated from the difference in 30-min NEE and R eco data (Eq. 7).
Gap-filled flux data are available through Ameri-Flux (http://ameriflux.lbl.gov) and additional processing details for eddy covariance data can be found in (Jimenez and others 2012;others 2014, 2016).

Statistical Analyses
We first created a matrix of all pairwise Pearson's correlation coefficients ( Figure S1 and S2) and examined univariate density plots and loess curves using the ggpairs function in the GGally package in R (Schloerke and others 2018) for the daily variables: NAPP (g C m -2 day -1 ), GAPP (g C m -2 day -1 ), R Aq (g C m -2 day -1 ), NEE (g C m -2 day -1 ), GEE (g C m -2 day -1 ), R eco (g C m -2 day -1 ), water level (m), air temperature (T air ;°C), and water temperature (T water ;°C) (Figures S1 and S2). We then evaluated the effect of water level and T water on daily NAPP, GAPP, and R Aq using generalized additive models (GAM) with the gam function in the ''mgcv'' package (Wood 2011) in R (R Core Team 2014). The GAM uses smooth functions to model nonlinear functional relationships between predictors and the response. Multiple predictors may be combined with the help of a tensor product smooth (Wood 2006(Wood , 2011. We included a tensor to smooth the interaction of water level and air temperature. This approach separates linear trends from any general nonlinear trends and determines if the significance of a smoothed variable is associated with a simple linear trend or a more complicated pattern (Wood 2011). We included water level, T water , and their interaction in the GAM analysis for all aquatic models and NAPP, water level, T air , and the interaction between water level and T air in ecosystem models. Temperature was included in models to capture the variation in thermal responses of the metabolisms, and water level was included to incorporate the size of the aquatic pool. No variables were removed from models to aid comparisons between ecosystems. We used Quantile-Quantile (Q-Q) plots of the residuals to evaluate if the appropriate distribution, identity and log-link were used. The Durbin-Watson statistic, with the durbinWatsonTest function in the ''car'' package (Fox and Weisberg 2011), was used to determine the AR(1) correlation coefficient. This coefficient was then used in the gam.check function to evaluate the convergence of the smoothness selection optimization and to run diagnostic tests to evaluate adequacy of dimension choices. We then proceeded by updating the original GAM model. The explanatory powers of the final models were compared using the R 2 statistic and the percent deviance explained. We also fit linear models to evaluate the relationship between the number of days inundated and annual ecosystem fluxes, NEE, GEE, and R eco . Because the assumption of independence between adjacent measurements was not met for annual data, this analysis is presented in a descriptive context for evaluation. All data used in this analysis is available on the Knowledge Network for Biocomplexity (Malone and others, 2021).

Ancillary Micrometeorological Measurements
Although seasonal fluctuations in water levels differed between the two sites, strong annual patterns in T air and PAR were consistent between the longhydroperiod marsh and the short-hydroperiod marl prairie over the study period (Figure 2). In the marsh, water levels remained above the soil surface throughout the study period. At the marl prairie, water levels dropped below the soil surface in the dry season of every year except in 2016. Over the study period, the mean water level in the marsh was 0.56 m and 0.13 m in the marl prairie. The mean number of dry days when water level was below the soil surface was 114 days per year in the marl prairie. The two sites were most similar during wet years (e.g., 2016; 0 dry days at the short-hydroperiod marl prairie) and differed during dry years (e.g., 2015; 248 dry days at the short-hydroperiod marl prairie).

Aquatic Ecosystem Metabolism
Temporal patterns in daily NAPP suggest that Everglades long-hydroperiod freshwater marshes were generally net heterotrophic, and mean daily absolute GAPP was 27% lower than the mean daily R Aq . Very few measurement days were net autotrophic (3 days). Daily GAPP ranged from 0 to -6.3 g C m -2 day -1 (mean = -0.7 g C m -2 day -1 ; Figure 3e), R Aq was between 0 to 6.1 g C m -2 day -1 (mean = 1.0 g C m -2 day -1 ; Figure 3i), and NAPP was between -0.1 to 0.8 g C m -2 day -1 (mean = 0.4 g C m -2 day -1 ). Patterns in NAPP were the result of shifts in GAPP and R Aq , which were differentially influenced by nonlinear relationships with water level, T water , and their interaction. Smoothed water level (p < 0.001), T water (p < 0.001), and their interaction (p < 0.001) explained variation in NAPP (R 2 = 0.31), GAPP (R 2 = 0.26), and R Aq (R 2 = 0.29) in the longhydroperiod marsh, suggesting there was a significant deviation from a flat line between explanatory variables and NAPP. The interaction between water level and T water showed that aquatic C uptake was greater than aquatic respiration when water levels were moderate ($ 0.5 m) and temperatures were high (Figure 3d). The short-hydroperiod marl prairie was consistently net heterotrophic ( Figure 4A) and mean daily absolute GAPP was 12% lower than R Aq . Daily GAPP ranged from 0 to -2.9 g C m -2 day -1 (mean = -0.7 g C m -2 day -1 ), R Aq was between 0 and 3.5 g C m -2 day -1 (mean = 0.8 g C m -2 day -1 ) and NAPP was between 0.1 and 0.7 g C m -2 day -1 (mean = 0.3 g C m -2 day -1 ). Water level (p < 0.001), T water (p < 0.001) and the interaction between the two (p < 0.001) explained variations in NAPP (R 2 = 0.33), GAPP (R 2 = 0.67), and R Aq (R 2 = 0.61) in the short-hydroperiod marl prairie, suggesting there was a significant deviation from a flat line between explanatory variables and NAPP. Similar to patterns observed in freshwater marsh, the interaction between water level and T water showed that aquatic C uptake was greater than aquatic respiration when water levels were moderate (0.2 to 0.4 m) and temperatures were high (Fig. 4d).

Ecosystem CO 2 Dynamics
At the ecosystem scale, the long-hydroperiod freshwater marsh was net heterotrophic (daily NEE > 0). Daily rates of NEE ranged from -1.1 to Figure 3. A, E, I Observed and generalized additive model predictions, and the effect of smoothed B, F, J water level, C, G, K water temperature (T water ) and D, H, I their interaction on A-D net aquatic primary productivity (NAPP), E-F gross primary productivity (GAPP), and I-L aquatic respiration (R Aq ) in the long-hydroperiod freshwater marsh (SRS-2). With atmospheric convention, C uptake is negative, and R 2 and the deviance explained (DEV) are used to evaluate the models.
1.2 g C m -2 day -1 (mean = 0.1 g C m -2 day -1 ), with most days (60%) exhibiting net heterotrophy (Fig. 5). Smoothed NAPP (p < 0.001), water level (p < 0.001), and the interaction between water level and T air (p < 0.001), explained variation in NEE (R 2 = 0.37), GEE (R 2 = 0.22), and R eco (R 2 = 0.64) in the long-hydroperiod marsh. Including NAPP in models increased the deviance explained by 5% for NEE, 4% for GEE, and 1% for R eco . NAPP had a positive relationship with NEE, suggesting that net heterotrophy in the water column coincided with net heterotrophy at the ecosystem scale at low rates of NAPP and that ecosystem heterotrophy was greatest at moderate NAPP (Figure 5b). The smoothed interaction between water level and T air showed that ecosystem fluxes were enhanced at moderate to high water levels (0.5-1.0 m) and at high and low temperatures (Figure 5e, j, o). Seasonal differences in C fluxes indicate that the GEE was greater than R eco during the wet season (Table 1). Annual NEE rates were variable in the long-hydroperiod marsh, with a mean annual value of 22 g C m -2 year -1 . Ecosystem respiration was 5% higher than GEE in the marsh, annually, and the number of inundated days had a strong correlation with annual R eco (R 2 = -0.33). The flux of C within the water column was within range of the ecosystem flux. Daily NAPP was 95% (± 151% standard errors) of daily NEE. Daily R Aq was 69% (± 2% standard errors) of R eco , and GAPP was 51% (± 2% standard errors) of daily GEE. It is important to note that the percentage of the ecosystem flux represented by the aquatic flux varied with water level without a clear directional relationship, and that data were also not equally representative across water level values.
At the ecosystem scale, short-hydroperiod marl prairies were also net heterotrophic (daily Figure 4. A, E, I Observed and generalized additive model predictions, and the effect of smoothed B, F, J water level, C, G, K water temperature (T water ) and D, H, I their interaction on A-D net aquatic primary productivity (NAPP), E-F gross primary productivity (GAPP), and I-L aquatic respiration (R Aq ) in the short-hydroperiod freshwater marl prairie (TS/Ph-1). With atmospheric convention, C uptake is negative, and R 2 and the deviance explained (DEV) are used to evaluate the models. Table 1. Total Seasonal Net Ecosystem Exchange (NEE; g C m -2 ), Gross Ecosystem Exchange (GEE; g C m -2 ), and Ecosystem Respiration (R eco ; g C m -2 ) at an Everglades Marsh (SRS-2) and Marl Prairie (TS/Ph-1) The wet season extends from May to October. With atmospheric convention, C uptake is negative NEE > 0) (Table 1). Most days (82%) exhibited net heterotrophy. Daily rates of NEE ranged from -3.1 to 3.2 g C m -2 day -1 (mean = 0.4 g C m -2 day -1 ) in the marl prairie (Figure 6a). Smoothed NAPP (p < 0.001), water level (p < 0.001), T water (p < 0.001) and the interaction between water level and T air (p < 0.001) explained variation in NEE (R 2 = 0.36), GEE (R 2 = 0.43), and R eco (R 2 = 0.21) in the short-hydroperiod marsh. Including NAPP in models increased the deviance explained by 4% for NEE, 1% for GEE, and 2% for R eco . Similar to patterns observed in the long-hydroperiod marsh, NAPP had a positive relationship with NEE, suggesting that net heterotrophy in the water column coincided with net heterotrophy at the ecosystem scale at low rates of NAPP but NEE leveled off at high rates of NAPP (Figure 6b). The smoothed interaction between water level and T air showed that ecosystem fluxes were enhanced at moderate to high water levels (0.2-4 m) and at high and low temperatures (Figure 6e, j, o). Seasonal differences in C fluxes indicated that GEE was greater than R eco during the dry season in the marl prairie. Annual NEE rates were variable, and the mean annual NEE was 133 g C m -2 year -1 in the marl prairie. Annually, R eco was 19% higher than GEE in the marl prairie and the number of inundated days per year had a strong correlation with annual R eco in the marl prairie (R 2 = 0.42). Daily NAPP was 64% (± 135% standard error) of daily NEE and 19% R eco (± 2% standard error). Daily R Aq was 56% (± 11% standard error) of R eco and GAPP was 53% (± 2 standard error) of daily GEE. Similar to patterns observed at the long-hydroperiod site, the percentage of the ecosystem flux represented by the aquatic flux was variable with water level without a clear directional relationship. The data was not equally distributed across water level values. With atmospheric convention, C uptake is negative, and R 2 and the deviance explained (DEV) are used to evaluate the models. Gaps in the modeled fluxes (A, F, K) are due to missing NAPP when the site was either dry or when data was missing.

DISCUSSION
We quantified the magnitude and variability of aquatic metabolic rates in a long-and a short-hydroperiod wetland and revealed their environmental controls. We also compared aquatic C dynamics to NEE and R eco in freshwater wetlands that differ in hydrology and the relative contribution of organic and inorganic soil C. Our results suggest that the complex interactions between water level and water temperature resulted in a high variability in aquatic metabolism that contributed to the net ecosystem CO 2 balance in subtropical freshwater marsh and marl prairie wetlands.
Inland aquatic ecosystems are primarily net heterotrophic (Duarte and Prairie 2005), and we measured higher respiration than photosynthesis within the water column throughout our study. In both the long-hydroperiod marsh and short-hydroperiod marl prairie wetlands, daily NAPP rates were similar to those of other freshwater ecosystems (Vaithiyanathan and Richardson 1998;Huggins 2005;Reeder 2011;Cooper and others 2013;Marois and others 2015;Kominoski and others 2018). In the water column, the productivity of algae and submerged aquatic vegetation are lower than the sum of heterotrophic processes that consume O 2 (Hagerthey and others 2010; Marois and others 2015), whereby more energy is being used than is being produced. Net heterotrophy is often observed in donor-controlled aquatic ecosystems that receive high amounts of allochthonous inputs of organic matter (Howarth and others 1996;Ram and others 2003;Cole and others 2006;Staehr and others 2012). Aquatic ecosystem respiration can be largely driven by substrate-specific microbial respiration of decomposing terrestrial-and aquaticderived organic matter (Webster and others 1999;Kominoski and others 2018).
We hypothesized that NAPP would increase with water level, indicating that aquatic metabolic activity would increase as the volume of the aquatic component of the system grew. We also hypothesized that aquatic productivity would increase with water temperature (Kemp and Testa 2011). Our results suggest that the interaction of water level and water temperature had an important effect on aquatic primary productivity. In both ecosystems, aquatic C uptake was greater than aquatic respiration simultaneously under moderate water levels, which represent the most frequent wet season condition, and at high temperatures. Patterns in NAPP indicate that above a water level of 0.5 m, increases in GAPP were smaller than R Aq at the long-hydroperiod freshwater marsh, perhaps due to changes in light attenuation. Although water levels in the Everglades increase in the summer wet season with PAR and water temperature contributing to an increase in the aquatic fluxes (Hagerthey and others 2010; Marois and others 2015), when floating periphyton is maximally developed, light penetration into the water column is greatly reduced (Online Appendix 1). In addition, the general enhancement of GAPP and R Aq with water level is likely due to the fact that wetlands with greater water depth and longer hydroperiods generally have higher aquatic species diversity (David 1996). Additionally, in the short-hydroperiod ecosystem, the longer the inundation persists, the greater the periphyton development is. In the higher-than-normal water levels of 2020, periphyton sweaters were 5-8 cm in diameter (Personal Observations).
We anticipated that NAPP would be positively correlated with NEE, indicating that net heterotrophy in the water column co-occurred and may have even facilitated the dominance of R eco . The significance of NAPP in models of NEE suggests that aquatic CO 2 dynamics in both long-and shorthydroperiod wetlands contribute to the vertical C flux. We also expected NAPP to have a stronger effect on the ecosystem fluxes at the long-hydroperiod marsh compared to the short-hydroperiod marl prairie because calcification in the water column may be weakening the relationship between aquatic and ecosystem fluxes in the shorthydroperiod marsh. Although the deviance explained by adding NAPP into the NEE model was very similar between sites, daily NAPP was a greater proportion of daily NEE at the long-hydroperiod site (95%) compared to the short-hydroperiod site (64%). Because the tower captures the slough 80% of the time in the long-hydroperiod site, the low representation of the emergent vegetation of the ridge is likely responsible for the higher contribution of NAPP to NEE. Although contributions to GEE were similar at long-and short-hydroperiod ecosystems, NAPP accounted for 45% of R eco at the long-hydroperiod site and 19% at the short-hydroperiod site.
In Everglades freshwater wetlands, there is a vertical separation of the net autotrophic and net heterotrophic components of the ecosystem. Emergent vegetation increases the autotrophic nature of wetland ecosystems and is a major source of the decoupling between NEE and NAPP. A chamber-based study in these Everglades long-and short-hydroperiod ecosystems found that emergent macrophytes dominated fluxes at the ecosys-tem scale (Schedlbauer and others 2012), whereas research focused on aquatic metabolism in this very landscape showed net heterotrophy within the water column (Hagerthey and others 2010). Shallow areas of aquatic ecosystems tend to be net autotrophic, often exporting particulate and dissolved organic matter to support respiration in adjacent regions of a water body (Kemp and others 1997;Caffrey and others 1998;Van de Bogert and others 2007). Past studies have shown how C produced in one region of an aquatic system may cause heterotrophy in adjacent regions following horizontal or vertical transport (Kemp and others 1999;Cole and others 2007;Williamson and others 2008;Lamberti and others 2010;Staehr and others 2012). Therefore, measuring only aquatic metabolism in ecosystems with emergent vegetation overestimates the heterotrophic nature of these ecosystems.
Understanding aquatic relative to ecosystem metabolism in wetlands is complex, especially in carbonate wetlands where C fluxes are associated with both metabolism and the carbonate cycle. Both NAPP and NEE may be differentially influenced by the process of calcium carbonate (CaCO 3 ) production and dissolution (Frankignoulle and others 1995;Macreadie and others 2017). Calcium carbonate production facilitates the release of CO 2 , which can either be used in organic C production through photosynthesis, or released to the atmosphere (Frankignoulle and others 1995;Zeebe and Wolf-Gladrow 2001;Pedersen and others 2013;Macreadie and others 2017). How aquatic vegetation uses CO 2 produced in the carbonate cycle for subsequent photosynthesis is not well-understood. However, the Everglades freshwater wetlands evaluated in this study can become a stronger source for CO 2 under inundation even when R eco should decline (Schedlbauer and others 2012;Zhao and others 2019), suggesting the carbonate cycle may be influencing the source potential of Everglades freshwater wetlands. CaCO 3 dynamics are a major factor in CO 2 fluxes in seagrass and carbonate-rich temperate coastal ecosystems (Barró n and others 2006). In Everglades ecosystems, the use of CO 2 from the carbonate cycle by aquatic vegetation (Allen and Spence 1981; Maberly and Spence 1983) may be causing a weakening of the relationship between aquatic metabolism and ecosystem metabolism which argues for more studies in this area. Continuous measurements of calcium carbonate production and dissolution in the longhydroperiod marsh compared to the short-hydroperiod marl prairie is needed to understand its impact on the decoupling of NAPP and NEE.
A substantive horizontal flux of C might add to the weak relationship we observed between aquatic and ecosystems fluxes. Studies have shown that estimates of lateral fluxes in upland to coastal ecosystems can be substantial (Krauss and others 2018). Direct and indirect estimates of exchange between wetlands (both marshes and mangrove forests) and coastal waters have demonstrated that lateral fluxes are often on the same order of magnitude as vertical exchanges others 2016, 2018;Ho and others 2017;Forbrich and others 2018;Krauss and others 2018;Maher and others 2018;Santos and others 2019). To date, studies have shown that the productivity of coastal wetlands serves a dual function of C burial and estuarine export, and the multiple fates of fixed C must be considered when evaluating wetland capacity for C sequestration (Krauss and others 2018). Given the complexity of the Everglades landscape, increasing the number of continuous, ecosystem-scale measurements of atmospheric CO 2 fluxes and a detailed accounting of biomass and soil elevation is needed to better understand the lateral flux of C across this landscape. Currently, comprehensive budgets including direct estimates of lateral C export are lacking (Bogard and others 2020).
Everglades wetlands slowly developed large organic C pools that are composed of two functionally different components, floc and peat (Gottlieb and others 2006;Craft and Richardson 2008;Pisani and others 2013). The highly labile and transient layer of floc formed from periphyton and SAV, accounts for the high rates of GAP and R Aq in both wetlands and contributes organic and inorganic materials to wetland soil and sediment formation (Pisani and others 2013). In contrast, peat consists mostly of refractory plant matter from above-and belowground compartments (cellulosic emergent and floating plants) that contributes to aerobic R Aq and R eco . Although studies suggest that some parts of the Everglades are sequestering C and accreting peat (Craft and Richardson 2008;Saunders and others 2015;Newman and others 2017), the vertical flux of CO 2 at our sites indicate that freshwater wetlands are currently net heterotrophic at the ecosystem scale annually others 2010, 2012;Jimenez and others 2012;Malone and others 2014). To better understand the sink/source capacity of Everglades freshwater ecosystems, we have to measure the horizontal flux of C.
Climate change simulations for subtropical ecosystems indicate that changes in the seasonal distribution of precipitation and small increases in annual precipitation will shift Everglades wetlands Integrating Aquatic Metabolism and Net Ecosystem to a more consistent C sink (Malone and others 2015), although this was based only on the vertical flux of C. This may be an indication that hydrological restoration is necessary to reduce the vertical loss of C, although more research is needed to evaluate the hydroperiods required to facilitate C sequestration and to understand the lateral flux of carbon across the landscape. As water level and temperature are related, water levels can have a strong moderating effect on respiration (Knox and others 2018). This relationship suggests that the warming induced vertical C loss could potentially be mitigated by water management practices. The restoration of nutrient-limited ecosystems like the Everglades may require decades for C storage capacity to increase. Legacies of hydrologic alteration, soil loss, phosphorus enrichment, and species invasions impact restoration efforts (the Comprehensive Everglades Restoration Plan) and ecosystem functioning (Newman and others 2017). For example, the temporal patterns in NEE that we measured are not representative of an ecosystem with high C storage capacity but this does not capture the lateral flux of carbon, limiting our understanding of the C budget.

Limitations to Data-Model Integration
Integration among ecosystem components remains a fundamental challenge in estimating wholeecosystem metabolism. In our study, breaks in time series due to seasonal dry downs and data loss limited our ability to evaluate lag effects and drivers of low-frequency changes in NAPP over time in both wetlands. Interpretation of the seasonal patterns of aquatic metabolism and relationships with environmental variables would be improved by information on changes in macrophyte and aquatic species composition and biomass. In particular, periphyton biomass is highly dynamic at both sites and is likely an important driver of net ecosystem CO 2 fluxes. In sloughs of the long-hydroperiod marsh, floating periphyton at peak development may cover most of the water surface, effectively shading deeper water. In marl prairies, periphyton requires several months to build up to peak biomass following the onset of the wet season, and at both sites periphyton degrades somewhat during the coolest months. In the sloughs of Shark River, Eleocharis and aquatic vegetation brown and die back with cooler weather. Measuring the links between seasonal changes in species activity and fluxes are necessary to understand the dynamic patterns of the integrated ecosystem metabolic balance in Everglades freshwater wetlands.
Integration of measurements and modeling across terrestrial-aquatic boundaries is required to understand and predict how global changes will continue to affect inland aquatic ecosystems. Largescale freshwater C process models lag behind aquatic nutrient (Wollheim and others 2008;Alexander and others 2009;Helton and others 2011) and terrestrial C models (for example, Terrestrial Ecosystem Model) (Melillo and others 1993), partly because of the complex nature of organic C pools, processing dynamics, and characterization of important drivers (for example flow, light, discharge) throughout networks. Understanding horizontal and vertical fluxes is also a major limiting factor in integrating aquatic dynamics in terrestrial C models. Moreover, a mismatch between the footprints of NEE and NAPP is a common issue for comparing variables measured with different approaches. For our study sites, the short-hydroperiod marl prairie generally had a homogeneous landscape while, in the longhydroperiod marsh the EC measurements were primarily from the slough (> 80% of the occasions according to the footprint analysis, data not shown) where NAPP was measured. Therefore, our conclusion may not be significantly influenced by this issue. For ecosystems with high heterogeneity, the footprint mismatch issue needs to be addressed by multiple NAPP measurements distributed throughout the footprint area of EC measurements.

CONCLUSIONS
We show critical links between aquatic metabolism and the same physical conditions that drive ecosystem CO 2 sequestration in wetlands, highlighting the importance of understanding NAPP in these systems where aquatic metabolism contributes to the C sequestering potential. These relationships are especially important in light of imminent changes in climate patterns, water management, and the Comprehensive Everglades Restoration Plan. The restoration of the natural hydrological patterns will be particularly important for controlling R eco , a driving force of the CO 2 sequestering capacity of Everglades freshwater wetlands. Ecosystems with high organic and inorganic C production are seldom examined to understand the role of these processes in ecosystem C sequestration (Barró n and others 2006). Studies often focus on organic C production by studying DO changes, which is only a good indicator of anaerobic respiration when the reduced products are oxidized (Barró n and others 2006). CaCO 3 production and dissolution prevents direct infer-ences on organic carbon fluxes from dissolved inorganic carbon (DIC) measurements. The examination of the metabolism and carbon fluxes requires the joint analysis of organic carbon and CaCO 3 fluxes, which is rarely measured simultaneously (Barró n and others 2006).
Understanding the integrated roles of C storage in aquatic and terrestrial ecosystems is becoming increasingly important, given the uncertainty in ability of coastal ecosystems to offset sea-level-rise with soil accretion. In many wetland ecosystems soil accretion rates are directly tied to aquatic productivity. The majority of inland aquatic ecosystems are net heterotrophic (Duarte and Prairie 2005), and global climate and land use changes are rapidly increasing net heterotrophy and C loss from ecosystems (Hotchkiss and others 2015; Rosemond and others 2015; Kominoski and others 2018). Quantifying the major sources and sinks of carbon in the biosphere to establish global carbon budgets requires a better understanding of carbon uptake, release and storage in aquatic ecosystems (Houghton 2007). Separating the aquatic from the ecosystem component allows us to better understand the relative impact of net autotrophic and net heterotrophic components of ecosystems, and differences in metabolic groups. The results of this study have implications for determining the influence of climate change and water management on the capacity to sequester C. Expanding our understanding of dynamic metabolic processes in wetland ecosystems across gradients in hydrology (Hagerthey and others 2010) will enhance our understanding of the role of aquatic networks in global elemental cycles (Benstead and Leigh 2012).

AC KNOWLEDGEMENT S
All research was performed under permits issued by Everglades National Park (EVER-2009-SCI-0070 and EVER-2013-SCI-0058). The authors would like to acknowledge the excellent support provided by the Florida Coastal Everglades Long-Term Ecological Research Program (FCE-LTER) and the Southeast Environmental Research Center in the Institute of the Environment at Florida International University. This manuscript is partly based upon work supported by the National Science Foundation (NSF) through the FCE-LTER program under Grant DEB-1237517. This research is also based in part on support from the Department of Energy's (DOE) National Institute for Climate Change Research (NICCR) through grant 07-SC-NICCR-1059 and the National Science Foundation Division of Atmospheric and Geospace Sciences Atmospheric Chemistry Program Awards 1561139, 1233006, 1801310 and 1807533. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of FCE-LTER, NSF, or DOE. The authors would also like to acknowledge the work of the reviewers who provided excellent advice on the data analysis and presentation of results, substantially improving the manuscript.

OPEN ACCESS
This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit h ttp://creativecommons.org/licenses/by/4.0/.

D A T A AV AI L ABI LI TY
Data can be found at https://knb.ecoinformatics.o rg/view/doi%3A10.5063%2FJQ0ZDM