The Unabated Atmospheric Carbon Losses in a Drowning Wetland Forest of North Carolina: A Point of No Return?

: Coastal wetlands provide the unique biogeochemical functions of storing a large fraction of the terrestrial carbon (C) pool and being among the most productive ecosystems in the world. However, coastal wetlands face numerous natural and anthropogenic disturbances that threaten their ecological integrity and C storage potential. To monitor the C balance of a coastal forested wetland, we established an eddy covariance ﬂux tower in a natural undrained bottomland hardwood forest in eastern North Carolina, USA. We examined the long-term trends (2009–2019) in gross primary productivity (GPP), ecosystem respiration (RE), and the net ecosystem C exchange (NEE) seasonally and inter-annually. We analyzed the response of C ﬂuxes and balance to climatic and hydrologic forcings and examined the possible effects of rising sea levels on the inland groundwater dynamics. Our results show that in 2009, a higher annual GPP (1922 g C m − 2 yr − 1 ) was observed than annual RE (1554 g C m − 2 yr − 1 ), resulting in a net C sink (NEE = − 368 g C m − 2 yr − 1 ). However, the annual C balance switched to a net C source in 2010 and onwards, varying from 87 g C m − 2 yr − 1 to 759 g C m − 2 yr − 1 . The multiple effects of air temperature (Tair), net radiation (Rn), groundwater table (GWT) depth, and precipitation ( p ) explained 66%, 71%, and 29% of the variation in GPP, RE, and NEE, respectively ( p < 0.0001). The lowering of GWT ( − 0.01 cm to − 14.26 cm) enhanced GPP and RE by 35% and 28%, respectively. We also observed a signiﬁcant positive correlation between mean sea level and GWT ( R 2 = 0.11), but not between GWT and p ( R 2 = 0.02). Cumulative ﬂuxes from 2009 to 2019 showed continuing C losses owing to a higher rate of increase of RE than GPP. This study contributes to carbon balance accounting to improve ecosystem models, relating C dynamics to temporal trends in under-represented coastal forested wetlands.


Introduction
Coastal wetlands serve important functions, such as providing food and bio-materials as direct resources, acting as wildlife sanctuaries, storing a large fraction of the terrestrial carbon (C) pool, protecting against storm surges, and facilitating sediment accumulation for land accretion [1][2][3]. They also provide water purification, tourism resorts, and other

The Study Site
The area of study is within the Alligator River National Wildlife Refuge (ARNWR), Dare County, North Carolina, USA, specifically at 35 • 45 16.44 N and 75 • 54 13.64 ( Figure 1). The eddy covariance flux tower is coded as US-NC4 in the FLUXNET registry. FLUXNET is a global network of micrometeorological tower sites that monitor the land-atmosphere exchange of carbon, water vapor, and energy fluxes. The US-NC4 flux tower was established in 2009. The tower stands within a >100-year-old natural coastal bottomland hardwood forest with natural drainage. Historically the site was used for timber production, and much of the surrounding area is farmland [54]. The refuge stretches across more than 154,000 acres of wetland. It is home to various wildlife (e.g., black bears, deer, etc.) and plant species (e.g., Nyssa sylvatica, Nyssa biflora, Taxodium distichum, Liquidambar styraciflua, and Pinus serotina). This wetland forest is also characterized by 'hummocks' around tree bases that are usually above the water table and non-vegetated low-lying 'hollows' submerged for more than 70% of the year [55]. The soil has 0.3-1.0 m thick organic soil horizons (Oi, Oe, and Oa), belonging to the soil type haplosaprist of the Pungo series (poorly drained with highly decomposed muck and less decomposed peat by highly reduced mineral sediments). The site lies <1 m above sea level. The average annual rainfall was 1163 mm , and the average temperature was 15.72 °C (2005-2017). More descriptions of the site can be found in [9,44,55,56].

Ecosystem Flux Measurement
An open-path infrared gas analyzer was used to measure the turbulent fluxes of CO2 exchange consisting of an LI-7500 (LI-COR, Lincoln, NE, USA), a Gill WindMaster sonic anemometer (R-350; Gill Instruments, Lymington, UK), and a CR1000 data logger, mounted at the top of the 32-m high tower (Fetch = ~2500 m).

Flux Data Processing and Gap-Filling
We calculated the net ecosystem exchange (NEE) for each 30-min period using the Eddypro software version 6.1.0 (LI-COR Inc., Lincoln, NE, USA). We used the standard quality checks and corrections for the following: spike detection [57], planar fit coordinate rotation of wind vectors [58], time lag corrections, and air density fluctuation [59], and The refuge stretches across more than 154,000 acres of wetland. It is home to various wildlife (e.g., black bears, deer, etc.) and plant species (e.g., Nyssa sylvatica, Nyssa biflora, Taxodium distichum, Liquidambar styraciflua, and Pinus serotina). This wetland forest is also characterized by 'hummocks' around tree bases that are usually above the water table and non-vegetated low-lying 'hollows' submerged for more than 70% of the year [55]. The soil has 0.3-1.0 m thick organic soil horizons (Oi, Oe, and Oa), belonging to the soil type haplosaprist of the Pungo series (poorly drained with highly decomposed muck and less decomposed peat by highly reduced mineral sediments). The site lies < 1 m above sea level. The average annual rainfall was 1163 mm (1981-2017), and the average temperature was 15.72 • C (2005-2017). More descriptions of the site can be found in [9,44,55,56].

Ecosystem Flux Measurement
An open-path infrared gas analyzer was used to measure the turbulent fluxes of CO 2 exchange consisting of an LI-7500 (LI-COR, Lincoln, NE, USA), a Gill WindMaster sonic anemometer (R-350; Gill Instruments, Lymington, UK), and a CR1000 data logger, mounted at the top of the 32-m high tower (Fetch =~2500 m).

Flux Data Processing and Gap-Filling
We calculated the net ecosystem exchange (NEE) for each 30-min period using the Eddypro software version 6.1.0 (LI-COR Inc., Lincoln, NE, USA). We used the standard quality checks and corrections for the following: spike detection [57], planar fit coordinate rotation of wind vectors [58], time lag corrections, and air density fluctuation [59], and high [60] and low pass filtering [61]. Flux outputs were flagged as 0 (high quality), 1 (medium), and 2 (low quality) when deviations in steady state and integral turbulence characteristics are <30%, >30% to <100%, and >100%, respectively [62].
During the post-processing of the 30-min data, we filtered the data for low signal strength, integral turbulence characteristics, despiking, and low friction velocity [63]. Missing NEE 30-min were gap-filled following the Max Planck Institute of Biogeochemistry method in Germany, which is an online protocol ReddyProc: Eddy covariance data processing tool (https://www.bgcjena.mpg.de/bgi/index.php/Services/REddyProcWeb) accessed on 20 February 2022. This web tool automatically partitioned the 30-min NEE into 30-min gross primary production (GPP) and ecosystem respiration (RE). The few remaining gaps (2%-6%) left after gap-filling and flux partitioning using ReddyProc [64] were filled by means of linear interpolation from the corresponding NEE values of the previous and succeeding years. Finally, daily NEE was computed as the sum of NEE 30min values over a 24-h period. Monthly averages were also obtained. The seasonality includes winter (January-March), spring (April-June), summer (July-September), and fall (October-December). The growing season usually occurs from April to October, while the non-growing season covers November to March. We present the annual C balance in atmospheric sign convention where positive NEE values connote C loss by the ecosystem, and negative values correspond to C gain [50].

Mortality and Leaf Area Index Determination
The mortality rates were determined from the 13 circular 7-m radius vegetation monitoring plots established within the flux tower footprint. Within each plot, tree species were identified and individual trees with diameter at breast height (i.e., 1.4 m aboveground level) greater than 2.5 cm were tagged. Diameter at breast height were measured each winter season. Dead trees were also recorded. The leaf area index (LAI) from 2009 to 2019 was generated using the Global Land Surface Satellite (GLASS) LAI datasets [65] from the Moderate Resolution Imaging Spectroradiometer (MODIS) surface reflectance data at an 8-day and 1-km resolution.

Data Analysis
We performed a generalized additive modeling using the gam function from the mgcv and MuMin packages to predict the best climate drivers of variations in C fluxes and balance [66]. In order to avoid over-parameterization, we used the second-order Akaike information criterion (AIC) in order to determine the best smoothing dimension with an upper limit equal to 3 [67]. The spline term effects were ranked according to their F-values (Wald Tests) [39]. The uncertainties in the C fluxes were obtained using the Monte Carlo bootstrapping approach. The Monte Carlo simulation was used to estimate a "probability" and the bootstrap approach was used to generate the same size of data from the same distribution many times (1000 times in this study). The bootstrap approach treated the original data as the population and was sampled with replacement from that data. The bootstrap function was used on this analysis [68]. Linear relationships were carried out using the ggplot2 [69], ggpubr, plotly [70], tidyverse [71], and reshape [72] packages. The confidence intervals of C fluxes were determined following [73] for the NEE and Monte Carlo approach using the bootstrap package used for GPP and RE. All analyses were processed in R version 4.1.1 [74].

Seasonal and Interannual Variations in Major Climate
Global radiation (Rg) and air temperature (Tair) followed a distinct seasonal trend of increasing after the winter season (March) or as the growing season started in spring (April), peaking in mid-summer (mid-July to mid-August) and gradually declining towards the non-growing season (November) ( Figure 2). However, groundwater table depth (GWT) and precipitation (p) had no apparent seasonality. The site was intermittently flooded in some periods of the year (mostly during non-growing seasons, max GWT = 24.18 cm), but frequently was also a little drier, mainly during the growing season (min GWT = −14.26)  depth (GWT) and precipitation (p) had no apparent seasonality. The site was intermittently flooded in some periods of the year (mostly during non-growing seasons, max GWT = 24.18 cm), but frequently was also a little drier, mainly during the growing season (min GWT = −14.26) (

Seasonal and Inter-Annual Variation in Carbon Fluxes and Balance
Across all years, the seasonal variation in GPP and RE were distinct. An increasing GPP and RE were apparent from the end of the non-growing season in winter towards spring. However, GPP and RE peaked in summer and gradually declined towards fall ( Figure 3). Average daily GPP and RE during the growing season were 8.09 g C m −2 d −1 and 8.31 g C m −2 d −1 , respectively. During the non-growing season, the average GPP was 2.25 g C m −2 d −1 , and RE was 3.04 g C m −2 d −1 . However, the seasonal variation in net ecosystem exchange was less distinct (Figure 3). During the growing season, NEE showed lower C loss (0.21 g C m −2 d −1 ) when the vegetation was active than during the non-growing season (0.79 g C m −2 d −1 ). Thus, photosynthesis and respiration during the growing season were higher than during the non-growing seasons, when GPP and RE were low ( Figure 3, p < 0.05). This pattern was fairly consistent across the study duration, except for anomalous spikes in 2010 and 2011, which we attribute to either flux system performance or random measurement errors.

Seasonal and Inter-Annual Variation in Carbon Fluxes and Balance
Across all years, the seasonal variation in GPP and RE were distinct. An increasing GPP and RE were apparent from the end of the non-growing season in winter towards spring. However, GPP and RE peaked in summer and gradually declined towards fall ( Figure 3). Average daily GPP and RE during the growing season were 8.09 g C m −2 d −1 and 8.31 g C m −2 d −1 , respectively. During the non-growing season, the average GPP was 2.25 g C m −2 d −1 , and RE was 3.04 g C m −2 d −1 . However, the seasonal variation in net ecosystem exchange was less distinct (Figure 3). During the growing season, NEE showed lower C loss (0.21 g C m −2 d −1 ) when the vegetation was active than during the non-growing season (0.79 g C m −2 d −1 ). Thus, photosynthesis and respiration during the growing season were higher than during the non-growing seasons, when GPP and RE were low ( Figure 3, p < 0.05). This pattern was fairly consistent across the study duration, except for anomalous spikes in 2010 and 2011, which we attribute to either flux system performance or random measurement errors.
In 2009, a higher annual GPP than the annual RE was observed (Table 1). This higher GPP than RE resulted in a net C sink. However, the annual carbon balance switched to a net C source from 2010 onwards, when RE continued to exceed GPP. The weak C source varied largely from 87 g C m −2 yr −1 in 2010 to as high as 759 g C m −2 yr −1 in 2017 (Table 1). We found a 27% inter-annual variation in GPP while that of RE was 36%. Forests 2022, 13, x FOR PEER REVIEW 6 of 19 In 2009, a higher annual GPP than the annual RE was observed (Table 1). This higher GPP than RE resulted in a net C sink. However, the annual carbon balance switched to a net C source from 2010 onwards, when RE continued to exceed GPP. The weak C source varied largely from 87 g C m −2 yr −1 in 2010 to as high as 759 g C m −2 yr −1 in 2017 (Table 1). We found a 27% inter-annual variation in GPP while that of RE was 36%.

Climatic Effects on Carbon Fluxes and Balance
The combined effects of Tair, Rn, GWT, and p explained 66%, 71%, and 29% of the variation in GPP, RE, and NEE, respectively (p < 0.0001, Table 2). Among the key climatic variables, the best fit models revealed Tair to be the main driver of GPP (F value = 119.29, p < 0.0001), RE (255.13, p < 0.0001), and NEE (42.31, p < 0.0001) ( Table 2). It was observed that water level fluctuations better explained GPP and RE than the radiation component, whereas p explained the least in the GPP, RE, and NEE models ( Table 2). Table 2. Results of generalized additive modeling (GAM) analysis to detect the best models and rank predictor variables e.g., air temperature (Tair), net radiation (Rn), groundwater table depth (GWT), and Precipitation (p) according to their ability to explain variation in gross primary productivity (GPP), ecosystem respiration (RE), and net ecosystem exchange (NEE). Daily data from 2009-2019 were used in the analysis. In order to avoid model over-parameterization, the second-order Akaike information criterion was used to determine the best smoothing dimension with an upper limit equal to 3. We report the coefficient of the fitted model with values for each spline term in turn. Since the variation in Tair and Rn remained relatively constant across years, we focus on GWT variation in the succeeding sections, as it was the most varied among the climate variables. It was also the second best predictor variable for GPP and RE ( Figure 2; Table 2). We evaluated how GWT was driving the carbon fluxes, especially RE, considering the prolonged inundation characteristic of the site, where tree mortality is accelerating.
When GWT levels were below the soil surface (GWT = −0.01 cm to −14.26 cm), GPP and RE were enhanced by 35% and 28%, respectively ( Figure 4). When GWT was taken as a sole predictor for GPP and RE alone, negative correlations were observed between GWT and GPP (R 2 = 0.19) and GWT and RE (R 2 = 0.20) during the growing season. However, the relationships between GWT with that of GPP and RE were weak during the non-growing season, with R 2 = 0.05 and R 2 = 0.02 for GPP and RE, respectively. (Figure 4). GWT had low explanatory power for the dynamics in NEE in both the growing and non-growing seasons.

Groundwater Level and the Rising Sea Level
GWT and mean sea level both increased over time, especially from 2015 onwards (Figure 5a), with a weak statistical relationship (R 2 = 0.11; Figure 5b). It is worth noting that the flux tower is approximately 20 km away from the coastline, which explains the weak relationship between GWT and mean sea level. We also observed a positive correlation between seawater level during tidal cycles and GWT (figure not shown). The relationship between GWT and sea level during low tide was higher (R 2 = 0.32) than during high tide (R 2 = 0.19).

Groundwater Level and the Rising Sea Level
GWT and mean sea level both increased over time, especially from 2015 onwards (Figure 5a), with a weak statistical relationship (R 2 = 0.11; Figure 5b). It is worth noting that the flux tower is approximately 20 km away from the coastline, which explains the weak relationship between GWT and mean sea level. We also observed a positive correlation between seawater level during tidal cycles and GWT (figure not shown). The relationship between GWT and sea level during low tide was higher (R 2 = 0.32) than during high tide (R 2 = 0.19).

The Cumulative Trends of Carbon Fluxes and Balance
A higher daily cumulative RE than GPP (p < 0.05) over the years of observation caused the ecosystem to be a net C emitter ( Figure 6).

The Cumulative Trends of Carbon Fluxes and Balance
A higher daily cumulative RE than GPP (p < 0.05) over the years of observation caused the ecosystem to be a net C emitter ( Figure 6).

The Cumulative Trends of Carbon Fluxes and Balance
A higher daily cumulative RE than GPP (p < 0.05) over the years of observation caused the ecosystem to be a net C emitter ( Figure 6).

Sources of Uncertainties
The eddy covariance measurement methodology is subject to uncertainties due to systematic errors [75] or random errors [76] that are difficult to assess. Despite the gapfilling methods developed for flux uncertainty, problems still exist [63,77,78]. Although these gap-filling errors affect the magnitude of the fluxes, all seasons and years should be affected similarly, allowing us to make comparisons over time.
In addition, our LAI is derived from remote-sensing data products. We caution the readers on using satellite-based LAI, as this approach may underestimate or overestimate LAI [79][80][81]. These uncertainties may occur due to overlapping or clumping between leaves and light obstruction from branches, boles, and stems [82].

Seasonal and Inter-Annual Variation in C Fluxes and Balance
The seasonal trends of NEE reflect the different phases and amplitudes in GPP and RE [83]. The largest carbon uptake occurred during the peak of the growing season (July-August) when irradiance was high, and temperature was at its peak. We observed a high GPP during the growing season, when greater solar radiation and high temperatures induced canopy development and high assimilation rates. This result is consistent with other studies [9,10,[84][85][86][87].
Nine years of monitoring the C balance of a natural coastal forested wetland along the southeastern US revealed how the ecosystem has transitioned from a net carbon sink during the measurement in 2009 into a net carbon source (2010-2017). Many old-growth forests have a decadal average of NEE close to zero [52] or are carbon sinks [88] unless the ecosystem is disturbed. Our site's long-term net C source condition may be ascribed to many tree physiological and environmental perturbations. These may include tree mortality, declining productivity, and thus LAI, prolonged hydroperiod, SLR, and saltwater intrusion, among others that are highly intertwined. We tackle each factor below with supporting data, except for saltwater intrusion and the construction of ditches, for which we based our explanations on previous work.

Effect of Tree Mortality on Carbon Fluxes
The decline in the annual rate of GPP compared to RE resulted in the ecosystem becoming a net carbon source over the years. A net loss of C to the atmosphere is likely due to the "drowning" of vegetation, causing increased tree mortality. Currently, there are hundreds of dead trees and stumps at ARNWR that were once part of a thriving wetland forest and are now part of a ghost forest. Ghost forests are made up of the remnants of the dead but still-standing trees. These trees are limbless and no longer have a leaf canopy or a functioning vascular system. We found only one dead tree per hectare a year after the start of measurements in 2009, yet that figure rose to 17 dead trees per hectare at the end of 2019 (Figure 7). Dead tree biomass may contribute up to 61% of all coarse woody debris on the forest floor [89], ultimately contributing to higher RE [90]. An increase of the number of dead trees on the forest floor indicates a major transformation of the ecosystem. This high incidence of mortality suggests that the environmental controls over NEE have slowly transitioned, resulting in the ecosystem losing stored C. A previous study from nearby the study site also reported 67% overstory mortality [27] and attributed this forest vegetation decline to increasing inundation. Over time, this increasing tree mortality leads to less canopy cover and shifts in species composition of the ecosystem.

Declining Leaf Area Index
LAI is a key driver of GPP [44]. In our study, inter-annual change in LAI was less variable, only ranging from 2.95 to 3.42. It must be noted that our wetland forest is a mature > 100-year-old ecosystem, and accordingly, LAI gradually stabilized with increasing age [91]. However, the significant reduction in LAI (R 2 = 0.70, Figure 7) from 2009-2019 suggests a decline in forest productivity over the years, thus turning the forest into a net carbon source.

Effect of Prolonged Hydroperiod on Carbon Fluxes
Relative to inundation, our site experienced a prolonged hydroperiod from 2009-2017, especially during the non-growing season (Figure 2). Extreme hydrological conditions alter wetland functioning [32]. High water level affects oxygen availability, soil organic matter decomposition rates, nutrient mobilization [85,92], and plant productivity [86,93]. During flooding, anaerobic respiration within the roots of plants produces toxic byproducts and limits the uptake of nutrients and water. Despite adaptations of some wetland trees to withstand inundation anoxia, constant flooding will cause stress since anaerobic conditions can only persist in the rooting zone for a few weeks during the growing season [94]. Flooding also triggers seedlings' mortality, limiting the population density, and thus, productivity. Also, plant litter decomposition generally responds quickly to changes in environmental conditions [95,96]. Litter decomposition may contribute most to respiration when the flooding subsides due to soil aeration and extracellular enzymatic ac-

Declining Leaf Area Index
LAI is a key driver of GPP [44]. In our study, inter-annual change in LAI was less variable, only ranging from 2.95 to 3.42. It must be noted that our wetland forest is a mature > 100-year-old ecosystem, and accordingly, LAI gradually stabilized with increasing age [91]. However, the significant reduction in LAI (R 2 = 0.70, Figure 7) from 2009-2019 suggests a decline in forest productivity over the years, thus turning the forest into a net carbon source.

Effect of Prolonged Hydroperiod on Carbon Fluxes
Relative to inundation, our site experienced a prolonged hydroperiod from 2009-2017, especially during the non-growing season (Figure 2). Extreme hydrological conditions alter wetland functioning [32]. High water level affects oxygen availability, soil organic matter decomposition rates, nutrient mobilization [85,92], and plant productivity [86,93]. During flooding, anaerobic respiration within the roots of plants produces toxic byproducts and limits the uptake of nutrients and water. Despite adaptations of some wetland trees to withstand inundation anoxia, constant flooding will cause stress since anaerobic conditions can only persist in the rooting zone for a few weeks during the growing season [94].
Flooding also triggers seedlings' mortality, limiting the population density, and thus, productivity. Also, plant litter decomposition generally responds quickly to changes in environmental conditions [95,96]. Litter decomposition may contribute most to respiration when the flooding subsides due to soil aeration and extracellular enzymatic activities of soil microorganisms. Prolonged flooding also affects seed germination and seedling recruitment [97]. Plants succumb to hydrological stress when it exceeds the tolerance limits of even flood-tolerant species. The continuing mortality at our site may lead to ecosystem transition from the mixed bottomland hardwood forest into a shrub-dominated pocosin with continued hydrologic stress, as is being seen in other parts of the Refuge [9,27].

Potential Impacts of Sea Level Rise
Projections reveal a future sea level ranging from +0.4 m to +1.2 m relative to today by 2100 [15,35]. Our analysis revealed that as SLR increases, the groundwater level around the study site is also increasing ( Figure 5). This increase of SLR indirectly decreased gross primary productivity and enhanced ecosystem respiration (Figure 4) through its effects on GWT, gradually shifting the forest towards a net positive carbon ecosystem exchange. This condition means that the forest will remain a net carbon source unless mitigated. Along with changes in carbon fluxes, we can expect to see a decrease of species richness and an overall decline in ecosystem diversity [21,28].
Rising sea levels might have impeded site drainage, lengthening the inundation period, and thus, tree mortality. This hydrologic stress might also be the reason for the small annual tree diameter growth and top die-back observed at our site. A previous record at the site from 2009-2019 for species such as bald cypress, red bay, water tupelo, and red maple showed a diameter increase of only 0.17 cm to 0.77 cm from 2009-2019 [9]. This phenomenon is not surprising since other wetland sites in Florida [98] and tidal wetlands in Louisiana [99] also observed this diameter growth suppression. Thus, continued monitoring is needed to examine the physiological problems of trees with prolonged inundation due to SLR.

Potential Impact of Saltwater Intrusion
In addition to hydrologic stress, saltwater intrusion may be a driving factor in the formation of ghost forests [28]. SLR may also increase saltwater intrusion into freshwater ecosystems [7,100], leading to species shifts and landward migration, and thus, to the loss of wetland areas [13,101]. Ghost forests result from exposure to the high salinity of water [28]. As saltwater seeps in, live trees first suffer from reduced annual growth. This manifests as a distressful condition, leading to death. As this continues, new seedling recruitment ceases [102,103]. Although we did not measure salinity exposure regimes in this study, the species at our site are vulnerable to salinity stress [97], and therefore, salinity needs to be monitored moving forward to disentangle the effects of salinity from hydrology. Casual spot measurements have never detected salinity at this relatively inland site, so we suspect the majority of the stress causing tree mortality and declining productivity is due to hydrology, but this is a hypothesis that is still in need of rigorous testing.

Major Climatic Effects on C Fluxes and Balance
The effect of solar irradiance, atmospheric air temperature, and soil water availability on GPP, RE, and NEE has been widely reported [39,40,87,[104][105][106]. It is important to determine the influence of these key climate drivers on forest C cycling because the climate-C cycle feedbacks control ecosystem-level C dynamics in response to the changing climate [45,52,53]. However, a deeper understanding of these climate effects on C balance remains limited, especially for coastal forested wetlands.
The combined effect of climatic drivers (Tair, Rn, GWT, and p) only explained 29% of the variability in NEE, suggesting greater biological control (e.g., LAI, mortality, etc.) on NEE than the climate drivers. Determining the response of carbon fluxes and balance to climate and vegetative factors is essential to understand the mechanisms behind the variations in GPP, RE, and NEE. However, obtaining biological data is more challenging than obtaining the climate variables [107]. Previous studies showed that biological effects could be larger [107,108], equivalent to [109], or less than [107,108] climatic effects at an interannual scale. Other studies also studied the correlations between climatic and biological effects [107,109], reflecting the responses of ecosystem C cycling to the changing climate.
Sustained favorable light availability resulting in increased GPP is well established in many studies [39,87,105,110] and is primarily a function of higher LAI [80,90]. However, at our site, GPP, RE, and NEE were not mainly controlled by net radiation but rather by temperature. This result suggests that, in this coastal wetland setting, caution must be taken in using solar irradiance to model C balance as other biological and environmental factors may contribute to the observed variability.
Air temperature is the major climatic control for RE for many studies [56,111,112]. Increased temperature enhanced litter decomposition in areas with abundant water [110,[113][114][115][116][117][118]. The sensitivity of RE to temperature may lead to continuing losses of C to the atmosphere due to hydrologic forcing coupled with a future warmer climate.
Surprisingly, GWT contributed more to the RE variability than Rn ( Table 2). The microtopography of our site influences the local hydrologic conditions. The hummocks surrounding tree bases are usually above the water table, while the hollows (usually nonvegetated) are submerged more than 70% of the year [55]. Site conditions such as this do not necessarily improve the GPP and RE. Instead, GPP and RE were enhanced during times of decreased groundwater level. Since NEE is the function of its component GPP and RE, the absence of variation in NEE with GWT suggests that even though GWT secondarily controls the variations in GPP and RE, they converge towards a relatively stable NEE. Therefore, modeling must be cognizant of using GWT to parameterize NEE in this mostly submerged wetland site.
During the non-growing season, the water table at our site is at full storage capacity [44], due to the low topographic setting and distinct microtopography. The site is a basin for the water from the nearby uplands. However, the water is difficult to drain because of the proximity of our site to the sea. The poor drainage characteristic of the site and poorly defined pathways for run-off decreased the overland flow. However, a deeper water table is attained during the growing seasons when precipitation is low. Thus, seasonal change in hydrology must be considered when modeling a unique coastal plain forest such as ours.

Why Do Changes in Groundwater Level Matter to Carbon Dynamics?
The site has low topographic setting characterized by seasonally dynamic high water levels. The variance explained for the negative linear relationship between GWT and GPP (R 2 = 0.19) and RE (R 2 = 0.20) appeared to behave similarly. However, a higher RE level (max = 29.78 g C m −2 d −1 ) than GPP (max = 27.45 g C m −2 d −1 ) (p > 0.05) when the GWT levels were decreased (GWT = 0 and below) suggests that the change in hydrology is affecting respiratory processes more than photosynthesis.
A previous study at our site reported that during non-flooded periods, 57% of the total ecosystem RE came from heterotrophic respiration, while 43% came from autotrophic respiration [55]. However, the reverse occurred during flooded periods. Autotrophic respiration contributed 69% of the total RE, while 31% came from heterotrophic respiration [55]. With this, we deduced that the greater sensitivity of ecosystem respiration to groundwater table drawdown was due to the higher sensitivity of soil respiration compared to plant respiration. Therefore, any reduction in water level at our site would positively stimulate soil respiration. More investigation is needed to better understand the hydrologic controls over RE in forested wetlands. Nevertheless, our results imply that increased hydrologic forcing and the changing climate could potentially shift the C balance trajectory of regional carbon fluxes and balance.

Is the Unabated Net C Source Strength in the Wetland Forest a Point of No Return?
Given continuing C losses for a prolonged period of time now (Figure 6), we expect that the formation of "ghost forests" brought about by the constant submergence of lowlying land will continue. Sea levels have risen at an average of 2.1-2.4 mm year −1 during the past two millennia [15,34,117]. Models suggest that coastal wetlands will migrate inland towards the freshwater wetlands, thereby increasing the wetland area up to 60% by the year 2100 for a 1.1 m SLR [21]. Carbon dating of organic sediments at a nearby site shows a historical organic soil accumulation rate of 1.11 to 1.13 mm year −1 , which is not enough to keep pace with the recent rise in sea level of 2.1 to 2.4 mm year −1 , resulting in ecosystem transition and ghost forest formation [26].
Therefore, a quantitative understanding of the effects of climatic extremes on coastal wetland hydrologic function is a pressing research need to inform wetland forest management. Studies that will determine the causes and consequences of forested coastal wetland ecosystem transitions are necessary to improve our understanding of the future of coastal wetland responses to environmental changes and the estimation of regional terrestrial C fluxes. Leaving this scenario as it is, without intervention, will finally set this ecosystem into a point of no return as the continental margin adjusts to SLR.

Forest Management Options
The accelerated transformation of wetland areas into ghost forests is alarming. In ARNWR alone, 32% of the refuge (31,600 hectares) has seen a change in land cover over 35 years, where 1151 ha of land is now open ocean, and about 19,200 ha of forest has transitioned to marshland or shrubland [107]. As of 2021, about 11% has transitioned to ghost forest [118].
Thus, wetland management interventions are needed to prevent aboveground carbon loss, whether through land preservation, reforestation, or the introduction of new tolerant species into the wetland environments. However, controversy arises in allowing new salt-tolerant vegetation to thrive and adapt within the wetland forests, as it may fail to preserve the current ecosystem.
What can be done to solve this issue and prevent above-ground carbon loss? While it may be too late to prevent harmful SLR, there are some mitigation strategies that should be used to limit the damage to this important coastal system. The most realistic strategy is to facilitate wetland migration. Essentially, this idea involves ensuring that salt tolerant vegetation is introduced as the ecosystem transitions. This prevents the area from turning into open water, and allows for the continued service of carbon sequestration. Ultimately, the long-term result would be the replacement of the ghost forest with a salt marsh ecosystem. Other strategies that should be explored include freshwater leaching and the engineering of structures that protect against SLR. However, these strategies are likely to be extremely costly, and are far less realistic than facilitated wetland migration. There is also a need for more research on the emergence of ghost forests and coastal forest resilience, particularly in regard to saltwater intrusion and salinity tolerance in these ecosystems. If nothing is done to mitigate this transition, North Carolina's wetland forest could be left with very little vegetation, decimating the ecosystem.

Conclusions
Over 11 years of monitoring C fluxes and balance in coastal forested wetland showed that the site had become a net carbon emitter for an extended period of time. The net C source strength was ascribed to the spike in tree mortality, declining productivity, prolonged hydroperiodicity, SLR, and potentially to saltwater intrusion. This natural wetland forest is highly affected by the hydrologic and climatic forcings that could worsen the future net C source trajectory, transitioning this ecosystem to a point of no return as the coastal margin adjusts to unabated SLR. Thus, various wetland management interventions or facilitative wetland migration plans must be explored to minimize aboveground carbon loss from these important coastal wetland forests in North Carolina.