Temporal Variations of Submarine Groundwater Discharge into a Tide-Dominated Coastal Wetland (Gaomei Wetland, Western Taiwan) Indicated by Radon and Radium Isotopes

: Submarine groundwater discharge (SGD) is evidenced around Taiwan, but the seasonal / temporal changes of SGD have not been fully examined. Here, we report a time-series investigation of SGD into a tide-dominated coastal wetland, the Gaomei Wetland, located to the south of the Da-Chia River’s mouth, western Taiwan, by using environmental tracers ( 222 Rn, 224 Ra ex , 228 Ra, δ D, and δ 18 O). Our results showed that regardless of dry and wet seasons, the 222 Rn activities in coastal waters were high at low tide but low at high tide. It represents the continuous input of 222 Rn-enriched groundwater. However, the 224 Ra ex and 228 Ra activities showed seasonal changes with tide conditions. In the dry season, the 224 Ra ex and 228 Ra activities in coastal waters were low at low tide but high at high tide; whereas in the wet season, an opposite relation was observed with quite high 224 Ra ex and 228 Ra activities in the low-tide waters. Coupled with the lower δ D and δ 18 O values of coastal and pore waters in the dry season, in comparison to those in the wet season, it is suggested that these phenomena probably reﬂected a seasonal di ﬀ erence in the main SGD component with fresh SGD in the dry season, but saline ones in the wet season. Based on a 222 Rn mass balance model, the estimated SGD ﬂuxes into the Gaomei Wetland varied with tidal ﬂuctuations and ranged from 0.2 to 25 cm d − 1 and from 0.1 to 47 cm d − 1 for the dry and wet seasons, respectively. A slightly high SGD ﬂux occurring during the wet season at spring tide, implied a stronger tidal pumping coupled with a larger hydraulic gradient between land and sea. In this study, we demonstrated that the variation of SGD into the Gaomei Wetland is not only controlled by the seasonal changes of groundwater recharge, but also by the tidal pumping process.


Introduction
Although not as obvious as surface runoffs, the flow of groundwater across the seabed into the ocean, i.e., Submarine Groundwater Discharge (SGD), does occur in many coastal areas. The first SGD observation was reported by Sonrel in 1868 [1], who described a submarine freshwater spring occurring at an 8 km-offshore site in the Mediterranean Sea, which the sailors in ancient times used to recharge freshwater on boat. As the highly spatial and temporal variation made it difficult to detect and measure SGD, the impact of SGD was neglected until the first elaborated estimate of SGD was Over the past decades, naturally occurring radionuclides, radon ( 222 Rn, T 1/2 = 3.82 d) and radium ( 223 Ra, T 1/2 = 11.4 d; 224 Ra, T 1/2 = 3.66 d; 226 Ra, T 1/2 = 1620 y; and 228 Ra, T 1/2 = 5.75 y) isotopes, have been used as tracers for detecting and quantifying local to regional scale SGDs as well as for studying their driving forces. Groundwater is greatly enriched in 222 Rn as compared to surface water and seawater (typically 1000-fold or greater), as it is continuously produced by the α decay of 226 Ra in the aquifers and it prefers to escape from water to an open environment like air [30]. Coupled with being chemically conservative and having a comparable half-life to residence time of coastal water, 222 Rn is an ideal tool for studying dynamic systems [30,31].
Radium (Ra) isotopes are produced by the radioactive decay of thorium (Th) isotopes, which are ubiquitous in sediments and rocks, and thought to be ideal tracers for RSGD. This is due to some of their characteristics, like their tendency to attach to the surface of sediment particles/rocks in freshwater (i.e., groundwater/fluvial water) but being released to water via ion exchange, which is induced when seawater mixes with freshwater in the estuarine/coastal zone, or is triggered by saline water recirculating in the coastal aquifer [31][32][33]. Accordingly, brackish waters are enriched in Ra isotopes; and when being transported away from the source region, activities of Ra isotopes are controlled by radioactive decay and mixing of water masses due to their conservation behavior in seawater.
Oxygen and hydrogen isotopic compositions (δ 18 O and δD) can further distinguish freshwater from saline water, due to a clear isotopic distinction between them [20,25,26,34]. During the hydrological process, continuous isotope fractionation causes the variation in δ 18 O and δD values in natural waters. As a result, seawater is significantly enriched in 18 O and D; while meteoric water as well as its subsequent surface runoff water and groundwater, are depleted in 18 O and D. Although regional and seasonal variations of the isotopic signature in meteoric water occur due to the different evaporation-precipitation conditions, the δ 18 O and δD values in natural waters, including atmospheric precipitation (meteoric water), surface runoff water, groundwater, brackish water, and seawater, tend to be linearly corrected [34][35][36][37][38]. When combined with other geochemical tracers in the coastal area, δ 18 O and δD serve as conservative indicators of the mixing between offshore seawater and onshore meteoric groundwater and further provide an insight to the source of SGD [35][36][37][38]

Study Area
The Gaomei Wetland is a tidal-dominated coastal wetland with a total area of 7.01 × 10 6 m 2 (3.5 km long and 1.8 km wide), which is connected to the Da-Chia River's estuary in the north as well as to the outlet of the Qingshui Drainage in the south (Figure 1), from where it also receives fresh fluvial waters. The Qingshui Drainage is an artificial system for supplying the agricultural water and draining the additional rainwater; while the Da-Chia River is a small mountainous river, originating from the Central Mountain Range, at an elevation of over 3500 m with a main stream length of 124 km that results in a mean relief ratio of 1/60 [39]. The drainage basin of the Da-Chia River, covering an area of 1235 km 2 , receives over 2000 mm of annual rainfall with more than 75% occurring during May-September; while the annual evaporation reaches up to 1200 mm [40]. All of these result in a high annual water discharge of 1.28 × 10 9 m 3 [39] into the Taiwan Strait and an annual groundwater recharge of over 500 mm into the highly permeable aquifers, hosted by unconsolidated gravelly/sandy sediments of the late Pleistocene and Holocene age [40]. In response to the topographic gradient in water tables between the recharge zone in the highland and the discharge zone in the lowland, the groundwater flows westward into the coastal plain of Taichung where the Gaomei Wetland is located, with a hydraulic gradient of~0.007 [41]. The surficial aquifer around the Gaomei Wetland is unconfined and about 40 m thick [41], consisting of Holocene alluvial gravelly/sandy sediments, with a high hydraulic conductivity of 6.18 × 10 −3~8 .3 × 10 −4 m s −1 . The surficial sediments of the Gaomei Wetland also consist of very fine sand and silt originating from the Da-Chia River, resulting in a relatively high permeability (10 −7~1 0 −4 m s −1 ; [41]). Semidiurnal tidal fluctuations with the greatest tidal ranges of over 4 m occurring at spring tide generate not only the flushing of seawater into the Gaomei Wetland, but also a significant and temporal change on sea levels. All of these imply a possible outflow of underground water across the seabed (i.e., SGD) into the overlaying water, and its magnitude might vary with tidal fluctuations. This flow is considered to be a pathway for discharging new/recirculated nutrients and other chemical materials into this wetland.

Field Work
Two field campaigns were carried out at the Gaomei Wetland under the contrasting hydrological conditions-during the dry season (9~10 May 2014; before the monsoon rainfall started) and during the wet season (27~28 August 2014; encountered a spring tide). We conducted time-series samplings over 30 h (over 2 tidal cycles) at a small-creek site 500 m offshore (Figure 1), where is fed by the fresh fluvial water and continuously flushed by the tidal seawater exchange. Water samples were taken every two hours at low tide and every one hour at high tide, for analyses of radon ( 222 Rn) and radium isotopes ( 224 Raex and 228 Ra), and oxygen and hydrogen isotopic compositions (δ 18 O and δD). In-situ temperature and salinity of water were recorded by the HORIBA U-50 series multiparameter water quality meter. In addition, pore water samples from 35 cm below the seabed were collected by using a push-point piezometer connected to a peristaltic pump, at a flow rate of ~200 mL min −1 . A simple weather station (Davis 6152C Cabled Vantage Pro2 Weather Station) with a propellertype anemometer mounted at a height of 4 m above the water surface, was set up at an inland site for recording on-site wind speed in 5 min interval.
For characterizing the water sources discharging into our study area, groundwater from the artesian wells around the Gaomei Wetland, river water from the mouth of the Da-Chia River and the outlet of the Qingshui Drainage were sampled ( Figure 1) in 2014 (August, December) and 2015 (September). They were analyzed for radon ( 222 Rn) and radium isotopes ( 224 Raex and 228 Ra), and oxygen and hydrogen isotopic compositions (δ 18 O and δD). For the seaward end-members, offshore surface and near-bottom seawater samples were collected at a site 3-km offshore (Figure 1; water depth ~28 m), using multiple 20-L Go-Flo bottles mounted on a SBE911+ CTD rosette assembly on 19 March 2014. The concurrent temperature and salinity were also recorded. The surficial aquifer around the Gaomei Wetland is unconfined and about 40 m thick [41], consisting of Holocene alluvial gravelly/sandy sediments, with a high hydraulic conductivity of 6.18 × 10 −3~8 .3 × 10 −4 m s −1 . The surficial sediments of the Gaomei Wetland also consist of very fine sand and silt originating from the Da-Chia River, resulting in a relatively high permeability (10 −7~1 0 −4 m s −1 ; [41]). Semidiurnal tidal fluctuations with the greatest tidal ranges of over 4 m occurring at spring tide generate not only the flushing of seawater into the Gaomei Wetland, but also a significant and temporal change on sea levels. All of these imply a possible outflow of underground water across the seabed (i.e., SGD) into the overlaying water, and its magnitude might vary with tidal fluctuations. This flow is considered to be a pathway for discharging new/recirculated nutrients and other chemical materials into this wetland.

Field Work
Two field campaigns were carried out at the Gaomei Wetland under the contrasting hydrological conditions-during the dry season (9~10 May 2014; before the monsoon rainfall started) and during the wet season (27~28 August 2014; encountered a spring tide). We conducted time-series samplings over 30 h (over 2 tidal cycles) at a small-creek site 500 m offshore (Figure 1), where is fed by the fresh fluvial water and continuously flushed by the tidal seawater exchange. Water samples were taken every two hours at low tide and every one hour at high tide, for analyses of radon ( 222 Rn) and radium isotopes ( 224 Ra ex and 228 Ra), and oxygen and hydrogen isotopic compositions (δ 18 O and δD). In-situ temperature and salinity of water were recorded by the HORIBA U-50 series multi-parameter water quality meter. In addition, pore water samples from 35 cm below the seabed were collected by using a push-point piezometer connected to a peristaltic pump, at a flow rate of~200 mL min −1 . A simple weather station (Davis 6152C Cabled Vantage Pro2 Weather Station) with a propeller-type anemometer mounted at a height of 4 m above the water surface, was set up at an inland site for recording on-site wind speed in 5 min interval.
For characterizing the water sources discharging into our study area, groundwater from the artesian wells around the Gaomei Wetland, river water from the mouth of the Da-Chia River and the outlet of the Qingshui Drainage were sampled ( Figure 1) in 2014 (August, December) and 2015 (September). They were analyzed for radon ( 222 Rn) and radium isotopes ( 224 Ra ex and 228 Ra), and oxygen and hydrogen isotopic compositions (δ 18 O and δD). For the seaward end-members, offshore surface and near-bottom seawater samples were collected at a site 3-km offshore ( Figure 1; water depth~28 m), using multiple 20-L Go-Flo bottles mounted on a SBE911+ CTD rosette assembly on 19 March 2014. The concurrent temperature and salinity were also recorded.

Measurement of 222 Rn in Water
Measurement of 222 Rn in water was carried out via using a RAD7-H 2 O system (Durridge Company, Inc., Billerica, MA, USA). Water samples were carefully taken into 250 mL glass bottles with no air bubbles occurring and sealed upon collection. After allowing 15 min for 218 Po to equilibrate with its parent 222 Rn, the sample bottle was connected to the RAD7-H 2 O system, equipped with an air-grabbing device. The 222 Rn air, together with 218 Po, was pumped into a closed air loop and detected through a radon-in-air detector. The RAD7 detector was calibrated using a standard prepared from a 226 Ra solution (Isotopes Products Laboratories, CF, Series No. 459-21-4), following the method of Xu et al. [42]. The activity of 222 Rn in water with 20% uncertainty was observable after a 30-min analysis, which was then decay-corrected to that of the sampling date and time. The limit of detection (LOD) for the RAD7-H 2 O system is 4 Bq m −3 .

Extraction and Measurement of Ra Isotopes in Water
For the determination of Ra isotopes, water samples (volume~40 L) were passed through a column packed with~25g manganese-coated acrylic fiber (MnO 2 fiber), at a flow rate of~1 L min −1 , to adsorb Ra isotopes onto the MnO 2 fiber. Based on the previous studies, over 95% of dissolved Ra could be retained on the fiber under the prescribed conditions [43,44]. Upon collection, each Ra-loaded fiber was rinsed with Ra-free water to remove sea salts, partially dried with a stream of compressed air, and sealed in a plastic bag for storage prior to laboratory procedures.
Within 3 days of collection, the Ra-loaded MnO 2 fiber was placed in a Radium Delayed Coincidence Counter (RaDeCC) system to count the α decays from 219 Rn and 220 Rn, daughter nuclides of 223 Ra and 224 Ra, respectively [45,46]. Afterwards the samples were recounted after being aged for 3~5 weeks, to allow the excess 224 Ra to decay and the supported 224 Ra to equilibrate with its parent 228 Th adsorbed onto MnO 2 fiber. In addition, the initial activity of the 228 Ra absorbed on the MnO 2 fiber was determined by measuring the sample, after being aged a least for 6 months to produce distinguishable 228 Th from 228 Ra decay [44]. The 224 Ra efficiency of the detector was determined by using a reference source prepared from thorium ores (IEAE-RGTh-1), in which the Ra isotopes were in secular equilibrium with their progenitors. The activities of 224 Ra and 228 Ra in water were calculated using the methods of Garcia-Solsona et al. [47] and Moore [43], with ±10% uncertainty.

Analyses of δ 18 O and δD of Water
Water sample was passed through a 0.45 µm filter into a 2-mL vial. This vial was immediately sealed with no existing air bubbles and was maintained at a room temperature, prior to analysis. Oxygen and hydrogen isotopic compositions (δ 18 O and δD) were analyzed by a cavity ring-down spectrometer (CRDS, Picarro L2140-i) housed at the Exploration and Development Research Institute, CPC Corporation, Taiwan. The results are reported as per mill (% ) relative to the Vienna Standard Mean Ocean Water (VSMOW). Replicate analyses δ 18 O and δD of the samples and variety of Standard Reference Materials (USGS- 45,46,47,48,49,50,53) showed an averaged 1 σ uncertainty better than 0.05% and 0.5% , respectively.

Station
For fluvial waters collected from the outlet of the Qingshui Drainage (salinity: 0.2~3.3) and the estuary of the Da-Chia River (salinity: 0~0.1), activities of 222 Rn were 554 ± 208 Bq m −3 and 1186 ± 305 Bq m −3 on average, respectively. These were much higher than those obtained from other rivers around Taiwan, such as the main stream of the Gaoping River in the southwest (264 ± 231 Bq m −3 ; [48]) and the Tanshui River in the north (271 ± 375 Bq m −3 ; [49]). The contribution of 222 Rn-enriched groundwater to fluvial water was suspected to be responsible for this, and this might also imply a significant exchange between groundwater and surface runoff water for the Da-Chia River. The fluvial water samples from the estuary of the Da-Chia River showed similar ranges of 224 Ra ex and 228 Ra activities to those of groundwater (Table 1); whereas in the outlet of the Qingshui Drainage, 224 Ra ex and 228 Ra activities were slightly high ( Table 1) and −50.0 ± 0.3% , respectively). This was attributed to a slight contribution of 18 O-and D-enriched offshore seawater.

Time-Series Observations of Coastal Water
Time-series observations on the coastal water at the small-creek site 500 m offshore, during the two field campaigns (9~10 May 2014 and 27~28 August 2014) are shown in Figure 2, together with the records of tidal height obtained from the tidal observatory (120 • 31 59" E, 24 • 17 16" N) [50] in the Port of Taichung, which is located to the south of the Gaomei Wetland.
Water 2020, 12, x FOR PEER REVIEW 7 of 22 Time-series observations on the coastal water at the small-creek site 500 m offshore, during the two field campaigns (9~10 May 2014 and 27~28 August 2014) are shown in Figure 2, together with the records of tidal height obtained from the tidal observatory (120°31′59′′E, 24°17′16′′N) [50] in the Port of Taichung, which is located to the south of the Gaomei Wetland.

9~10 May 2014 (Dry Season).
As this field campaign was conducted before the rainy season started, we considered these results to be representative for the dry season. In addition, we encountered a neap tide, thus, the tidal range during our sampling period was much less than 3 m. This implied a relatively low contribution of the offshore seawater to the Gaomei Wetland during high tide, resulting in a slightly lower salinity of surface water (salinity: 24.6, consisting of ~70% seawater); while during low tide, freshwater (salinity: ~0.35) occupied the tidal flat. The water

Variations of Hydrographic Features
9~10 May 2014 (Dry Season). As this field campaign was conducted before the rainy season started, we considered these results to be representative for the dry season. In addition, we encountered a neap tide, thus, the tidal range during our sampling period was much less than 3 m. This implied a relatively low contribution of the offshore seawater to the Gaomei Wetland during high tide, resulting in a slightly lower salinity of surface water (salinity: 24.6, consisting of~70% seawater); while during low tide, freshwater (salinity:~0.35) occupied the tidal flat. The water temperature fell into a range of 21.8~29.0 • C. The wind speed ranged between 0.0~4.0 m s −1 , with high ones occurring in the afternoon.
27~28 August 2014 (Wet Season). Since we encountered a spring tide during this field campaign, the tidal range was more than 4 m, which meant that more seawater intruded into the Gaomei Wetland during the flood tide. As a result, seawater (salinity:~34.5) covered the tidal flat at high tide and a part of them remained during the low-tide period, resulting in an average salinity of~22.0 (consisting of 63% seawater). The water temperature ranged between 27.0~33.7 • C and the wind speed ranged between 0.0~3.6 m s −1 , with higher ones prevailing during night-time. 27~28 August 2014 (Wet Season). Activities of the 222 Rn in coastal waters varied from~0 to 548 Bq m 3 , with an average of 181 ± 173 Bq m −3 . These fluctuated with tidal cycles, displaying high activities during low-tide period but low ones during high-tide periods and reaching to peak activities about every 12 h. The activities of 224 Ra ex ranged from 3.33 to 49.3 dpm 100 L −1 and showed a contrasting trend to those during the dry season. Coastal waters with 224 Ra ex activities higher than 20 dpm 100 L −1 were found during the low-tide period with their salinity ranging between 19.1 and 23.4; while low 224 Ra ex activities (<10 dpm 100 L −1 ) occurred during the high-tide period. The temporal variation of 228 Ra activities also showed the same trend to those of 224 Ra ex and fell into a range of 3.68~65.4 dpm 100 L −1 . The δ 18 Oand δD values of coastal waters showed a similar temporal variation to those at dry season and ranged from −8.48 to −0.10% and from −59.4 to −1.1% , respectively.

Time-Series Observations of Pore Water
Time-series sampling of pore water from 35 cm below the seabed were also carried out at the small-creek site 500 m offshore in the Gaomei Wetland, during the two field campaigns (9~10 May 2014 and 27~28 August 2014), and was analyzed for salinity, 222 Rn, δ 18 O, and δD. The results are plotted in Figure 3 along with the tidal height.
9~10 May 2014 (Dry Season). Salinity of pore water varied from 3.70 to 7.35. Activities of 222 Rn in pore waters varied from~0 to 1526 Bq m 3 with an average of 558 ± 137 Bq m −3 , slightly higher than those in the coastal waters. These fluctuated with tidal cycles, but displayed different patterns from those in the coastal water. Interestingly, 222 Rn activities in pore waters increased during the flooding period but decreased during the ebbing period, reaching a peak activity in about 6 h after low tide. The δ 18 Oand δD values ranged from −6.25 to −5.16% and from −43.0 to −38.3% , respectively, which were just slightly higher than those in the groundwater samples.
27~28 August 2014 (Wet Season). Salinity of pore water varied from 10.67 to 22.62, which were higher than those during dry season, displaying a 3-h delay with tidal fluctuations. Activities of 222 Rn in pore water ranged from~0 to 1323 Bq m 3 , with an average of 566 ± 147 Bq m −3 , showing the same pattern to that of dry seasons. Although 222 Rn activities increased during the flooding period but decreased during the ebbing period, they reached a peak activity in about 3 to 4 h after the low tide. The lag time between low tide and the peak of the 222 Rn concentration in pore water was possibly due to the longer time required to fill the shallow aquifer with 222 Rn-riched groundwater than to drain it out. The δ 18 O and δD values varied between −3.98 and −0.37% and between −27.3 and −2.3% , respectively, was found to be relatively higher, in comparison to those at dry season, and showed the same trend in salinity. All of this indicates that seawater re-penetrating into the seabed is equilibrated in 3 h after the high tide.

Observations of δD and δ 18 O
Atmospheric precipitations are the main source of groundwater recharge in the Taichung groundwater catchment (i.e., the hinterland of the Gaomei Wetland). Thus, the isotopic signatures recorded in groundwater should be obedient to those in meteoric water. The Local Meteoric Water Line (LMWL), a least-squares regression line from the stable isotopic compositions in precipitation, as recommend by the International Atomic Energy Agency (IAEA) [51], is composed of long-term isotopic signatures of the local precipitation. The calculated LMWL for Taichung is represented as follows [52]: δD = 7.9 × δ 18 O+12.8 (1) which is parallel to the Global Meteoric Water Line (GMWL) defined by Craig [53]: In the δD versus δ 18 O diagram we plotted with all data determined in this study (Figure 4), the distributions of groundwater and fresh fluvial water samples fit well with the LMWL for the Taichung groundwater catchment. This indicates the following two features-(1) the stable isotopes in groundwater and fresh fluvial water were unaffected by water-rock interactions or significant evaporation and could be regarded as conservative; (2) coupled with a remarkably high 222 Rn activities (1,186 ± 216 Bq m −3 ) in fluvial water sample, we suggest there might be a contribution of groundwater to fluvial water. This further implied a significant exchange between underground water and surface runoff water for the Da-Chia River. In contrast, the isotopic compositions of brackish fluvial water, coastal water, pore water, and seawater distributed below the LMWL line, that showed a regression line (δD = 6.7 × δ 18 O−0.9; R 2 = 1.00) with a slightly lower slope, as compared to the LMWL. This distinguishing feature could be due to the seawater-atmospheric water vapor interactions [36]. It further indicates that the coastal and pore waters at the Gaomei Wetland originated from meteoric water and then mixed with offshore seawater.

Observations of δD and δ 18 O
Atmospheric precipitations are the main source of groundwater recharge in the Taichung groundwater catchment (i.e., the hinterland of the Gaomei Wetland). Thus, the isotopic signatures recorded in groundwater should be obedient to those in meteoric water. The Local Meteoric Water Line (LMWL), a least-squares regression line from the stable isotopic compositions in precipitation, as recommend by the International Atomic Energy Agency (IAEA) [51], is composed of long-term isotopic signatures of the local precipitation. The calculated LMWL for Taichung is represented as follows [52]: which is parallel to the Global Meteoric Water Line (GMWL) defined by Craig [53]: In the δD versus δ 18 O diagram we plotted with all data determined in this study (Figure 4), the distributions of groundwater and fresh fluvial water samples fit well with the LMWL for the Taichung groundwater catchment. This indicates the following two features-(1) the stable isotopes in groundwater and fresh fluvial water were unaffected by water-rock interactions or significant evaporation and could be regarded as conservative; (2) coupled with a remarkably high 222 Rn activities (1186 ± 216 Bq m −3 ) in fluvial water sample, we suggest there might be a contribution of groundwater to fluvial water. This further implied a significant exchange between underground water and surface runoff water for the Da-Chia River. In contrast, the isotopic compositions of brackish fluvial water, coastal water, pore water, and seawater distributed below the LMWL line, that showed a regression line (δD = 6.7 × δ 18 O − 0.9; R 2 = 1.00) with a slightly lower slope, as compared to the LMWL. This distinguishing feature could be due to the seawater-atmospheric water vapor interactions [36]. It further indicates that the coastal and pore waters at the Gaomei Wetland originated from meteoric water and then mixed with offshore seawater.  Figure 5A), which was similar to that observed off the coast of Southern Taiwan (δ 18 O = 0.24 × S − 8.33, R 2 = 0.97; [26]). This confirmed that groundwater, even pore water samples, were well separated from seawater. A consistent pattern was shown in the plot of δD versus salinity (δ 18 O = 1.72 × S − 57.36; R 2 = 0.88; Figure  5B). These further indicated that the variability in the δ 18 O and δD values was induced by the different fraction of seawater contributing to the collected coastal water and pore waters. The fraction of seawater in pore water could be simply calculated via two end-member mixing model of salinity. We observed about 11%~65% of seawater in the collected pore water (35 cm below the seabed), implying that the re-circulated seawater played an important role in groundwater-seawater interactions in the coastal aquifer, which meant that the re-circulated seawater mixed with fresh groundwater (i.e., brackish groundwater) might be the dominated source for SGD in the Gaomei Wetland.   Figure 5A), which was similar to that observed off the coast of Southern Taiwan (δ 18 O = 0.24 × S − 8.33, R 2 = 0.97; [26]). This confirmed that groundwater, even pore water samples, were well separated from seawater. A consistent pattern was shown in the plot of δD versus salinity (δ 18 O = 1.72 × S − 57.36; R 2 = 0.88; Figure 5B). These further indicated that the variability in the δ 18 O and δD values was induced by the different fraction of seawater contributing to the collected coastal water and pore waters. The fraction of seawater in pore water could be simply calculated via two end-member mixing model of salinity. We observed about 11%~65% of seawater in the collected pore water (35 cm below the seabed), implying that the re-circulated seawater played an important role in groundwater-seawater interactions in the coastal aquifer, which meant that the re-circulated seawater mixed with fresh groundwater (i.e., brackish groundwater) might be the dominated source for SGD in the Gaomei Wetland. Water 2020, 12, x FOR PEER REVIEW 11 of 22  Figure 6 shows the diagrams of salinity versus 224 Raex and 228 Ra with a gray dash line representing the expected activity due to the conservative mixing between offshore seawater, and fresh fluvial water from the Da-Chia River and the Qingshui drainage. Except for the coastal water samples collected during the low-tide period of the dry season, all data points were distributed above the expected activity line, especially for the samples collected during the low-tide period of the wet season (salinity: 19.1~23.4), which were significantly higher than the expected. The excess concentrations of 224 Raex and 228 Ra would have originated from other sources, such as SGD, desorption the riverine suspended particles, and diffusion from benthic sediments. Here, quite low contributions (<5%) via diffusion from the benthic sediments to the total 224 Ra and 228 Ra budget were typically reported, thus its contribution could be ignored [54,55].   Figure 6 shows the diagrams of salinity versus 224 Ra ex and 228 Ra with a gray dash line representing the expected activity due to the conservative mixing between offshore seawater, and fresh fluvial water from the Da-Chia River and the Qingshui drainage. Except for the coastal water samples collected during the low-tide period of the dry season, all data points were distributed above the expected activity line, especially for the samples collected during the low-tide period of the wet season (salinity: 19.1~23.4), which were significantly higher than the expected. The excess concentrations of 224 Ra ex and 228 Ra would have originated from other sources, such as SGD, desorption the riverine suspended particles, and diffusion from benthic sediments. Here, quite low contributions (<5%) via diffusion from the benthic sediments to the total 224 Ra and 228 Ra budget were typically reported, thus its contribution could be ignored [54,55].  Figure 6 shows the diagrams of salinity versus 224 Raex and 228 Ra with a gray dash line representing the expected activity due to the conservative mixing between offshore seawater, and fresh fluvial water from the Da-Chia River and the Qingshui drainage. Except for the coastal water samples collected during the low-tide period of the dry season, all data points were distributed above the expected activity line, especially for the samples collected during the low-tide period of the wet season (salinity: 19.1~23.4), which were significantly higher than the expected. The excess concentrations of 224 Raex and 228 Ra would have originated from other sources, such as SGD, desorption the riverine suspended particles, and diffusion from benthic sediments. Here, quite low contributions (<5%) via diffusion from the benthic sediments to the total 224 Ra and 228 Ra budget were typically reported, thus its contribution could be ignored [54,55].  In order to calculate the excess Ra concentrations in the coastal water, the chemical mass balance for 228 Ra and 224 Ra ex could be expressed by the following equations:

Excess 224 Raex and 228 Ra in the Coastal Water
where Ra M , Ra S , and Ra R are the Ra activities in the coastal water, offshore seawater ( 224 Ra ex : 1.96 ± 0.33 dpm 100L −1 , 228 Ra: 1.34 ± 0.21 dpm 100L −1 ), and fresh fluvial water ( 224 Ra ex : 2.55 ± 0.47 dpm 100L −1 , 228 Ra: 2.85 ± 0.37 dpm 100L −1 ) end-members, respectively; f S is the fraction of offshore seawater, which was calculated from the ratio of the measured salinity to that of the offshore water (salinity: 34.7). According to the reported desorption experiment of Ra [56,57], the ion-exchange equilibrium was achieved within 0.5~4 h after being exposed to saline water and it reached the maximum desorption of Ra in the salinity range of 10~20. Furthermore, about 22~41% of Ra on riverine suspended sediment were available for ion-exchange [58]. In Table 2, the maximum dissolved 224 Ra and 228 Ra activities via desorption from riverine suspended particles were roughly calculated and the results were 12.1 ± 2.13 and 12.5 ± 2.22 dpm 100L −1 , respectively. Subsequently, the excess 224 Ra ex and 228 Ra contributed by desorption from riverine suspended particles to coastal water could be subtracted by the following equation [58]: where (E xcess 228 Ra) and (E xcess 224 Ra ex ) are the corrected excess 228 Ra and 224 Ra ex activities in the coastal water sample; Ra DES is the maximum Ra activity desorbing from riverine suspended particles (Table 2). Sal S is the salinity of coastal water sample and Sal U is the salinity at which the ion-exchange equilibrium is achieved. Here, we took the salinity of 20 for Sal U . During the dry season, the corrected excess 228 Ra and 224 Ra ex activities in coastal water were insufficient during the low-tide period, but fell into 0.52~4.42 and 1.98~10.13 dpm 100 L −1 , respectively, during the high-tide period. During the wet season, they ranged from 1.24 to 60.4 dpm 100 L −1 and from 0.17 to 44.4 dpm 100 L −1 for 228 Ra and 224 Ra ex , respectively, with significantly high values in low-tide waters. We suspected that such extremely high excess 224 Ra ex and 228 Ra activities in low-tide waters during the wet season might be caused by other 224 Ra-and 228 Ra-enriched inputs. Although measurements of 224 Ra and 228 Ra in pore water were not carried out in this study, we inferred that RSGD might be a potential source, based on the signatures of δD and δ 18 O in pore waters. 0.71 0.74 Riverine suspended particle (g 100L −1 ) [39] 17 ± 3 17 ± 3 Maximum dissolved Ra by desorption from suspended particle (dpm 100L −1 ) ** 12.1 ± 2.13 12.5 ± 2.22 * Total 224 Ra in sediments was not measured directly, but was inferred from the 228 Th. ** The calculation followed Krest et al. [58].

222 Rn mass Balance for SGD Rate Estimation
The change in the concentration of 222 Rn in coastal water of the Gaomei Wetland over time is a function of (1) production from 226 Ra dissolved in the seawater; (2) benthic input from sediments, including both SGD and diffusion; (3) input of 222 Rn from the Da-Chia River; (4) loss of 222 Rn via atmospheric evasion; (5) lowering of the 222 Rn concentration by mixing with lower concentration seawater offshore; and (6) loss of 222 Rn via decay. For quantifying SGD in the Gaomei Wetland via time-series 222 Rn observations, a mass balance (Figure 7) was constructed as follows [30,[60][61][62]: where F(t) is the net flux of 222 Rn at a specific time (t); F SGD , F Ra-226 , and F River represent the 222 Rn flux attributed to SGD, supported by 226 Ra in water, and discharged by the Da-Chia River, respectively. F Diff is the diffusion flux of 222 Rn from the bottom sediments; F atm is the 222 Rn flux to the atmosphere, and F decay is the decay loss of 222 Rn flux. F T-in and F T-out represent the tidal effect on 222 Rn inventory, which were the fluxes entering and leaving with the incoming and outgoing tides, respectively.
( ) where F(t) is the net flux of 222 Rn at a specific time (t); FSGD, FRa-226, and FRiver represent the 222 Rn flux attributed to SGD, supported by 226 Ra in water, and discharged by the Da-Chia River, respectively. FDiff is the diffusion flux of 222 Rn from the bottom sediments; Fatm is the 222 Rn flux to the atmosphere, and Fdecay is the decay loss of 222 Rn flux. FT-in and FT-out represent the tidal effect on 222 Rn inventory, which were the fluxes entering and leaving with the incoming and outgoing tides, respectively.
If all terms are known, the FSGD can be calculated and converted to SGD flux (QSGD) by dividing by the mean 222 Rn concentration in groundwater collected from artesian well, in the vicinity of the Gaomei Wetland (CGW: 9,187 ± 595 Bq m −3 ): The detailed calculation for each term in the 222 Rn mass balance model are described in the following subsections.  If all terms are known, the F SGD can be calculated and converted to SGD flux (Q SGD ) by dividing by the mean 222 Rn concentration in groundwater collected from artesian well, in the vicinity of the Gaomei Wetland (C GW : 9187 ± 595 Bq m −3 ): The detailed calculation for each term in the 222 Rn mass balance model are described in the following subsections.

The Temporal Variations of net 222 Rn fluxes
To determine the 222 Rn flux associated with SGD in the Gaomei Wetland, the net 222 Rn flux at each time (t) was first calculated through the following equation, under the assumption of well-mixing water column [61]: where T f represents the flushing time of water in the Gaomei Wetland (0.64 d; detail calculation will be described in the "Tide effect" subsection) and λ is the decay constant of 222 Rn (7.56 × 10 −3 h −1 ). I(t), the 222 Rn inventory at each time (t), is defined as the product of the measured concentration of 222 Rn ( 222 Rn M ) and water depth (h) at each time (t): It is regarded as 222 Rn in a coastal-water column within an area of 1 m 2 . Here, we subtracted the activity of 222 Ra in offshore water ( 226 Ra S : 0.83 Bq m −3 ) from the measured concentration of 222 Rn ( 222 Rn M ) to obtain the excess (unsupported) 222 Rn. The estimates of net 222 Rn flux at each time (t) were between 23.0 and 108 Bq m −2 h −1 for the dry season and between 14.0 and 192 Bq m −2 h −1 for the wet season.

Atmospheric loss
Gas exchange between air and water always leads to 222 Rn loss to the atmosphere, and the flux could be described as follows [63]: where F atm is the diffusion flux across the water-air interface (Bq m −2 h −1 ); α represents the Ostwald's solubility coefficient, which is temperature-dependent (α = 0.105 + 0.405e −0.0502T , T is the temperature of water in • C; [30]); C w and C a are the activity of 222 Rn in water and air, respectively. Here, C w for 222 Rn was a measure of the RAD7-H 2 O system; while we took 5 Bq m −3 as C a for 222 Rn in Taiwan [64].
In addition, the parameter, k, represents the gas transfer velocity (m s −1 ), which is the key parameter for the diffusion flux across the water-air interface. An empirical equation relating k to temperature and wind speed is proposed as follows [63]: k 600 = 0.45U 10 1.6 (Sc/600) −a (12) where U 10 is the wind speed at 10 m height above the water surface; a is a variable exponent that equals 0.6667 for U 10 ≤ 3.6 m s −1 and equals 0.5, when U 10 > 3.6 m s −1 . Our in-situ wind speeds, which were recorded at the height of 4 m above the water surface, were converted to those at the height of 10 m, via the following equation reported by Pond [65]: where z is the observation elevation (m) and a = 0.097. Moreover, the parameter, Sc, represents the Schmidt number, which is a ratio of the kinematic viscosity (ν) of water to the molecular diffusion coefficient (D m ) of gas in water, Sc = v/D m , and was divided by 600 to normalize k to CO 2 at 20 • C in freshwater. The molecular diffusion coefficient could be calculated via a temperature-dependent function [66,67]: where T is the temperature of water in • C. The kinematic viscosity (ν) is a simple ratio of the absolute viscosity (µ) to the density (ρ) of water at a measured temperature, v = µ/ρ. In-situ density of water was recorded by the HORIBA U-50 series multi-parameter water quality meter. To calculate the absolute viscosity, we employed an empirical relationship associated to water temperature [68]: where the coefficients in the equation: A = 2.414 × 10 −5 kg s −1 m −1 ; B = 247.8 K; C = 140 K; and T represents water temperature ( • C). Based on the aforementioned method, the diffusion fluxes of 222 Rn across the water-air interface ranged from 0.02 to 8.8 Bq m −2 h −1 during the wet season and from 0.01 to 14.3 Bq m −2 h −1 during dry seasons. Overall, the results showed slightly high atmospheric evasion during low tide and vice versa during high tide, because of higher concentration in water during low tide than that during high tide. This phenomenon indicated that the difference of 222 Rn concentration in water and air might a dominating factor in the flux of atmospheric evasion.

Diffusion Flux from Bottom Sediments
Rn-222 will diffuse from the bottom sediments into the overlying water because the 222 Rn concentrations in pore waters are much higher. The diffusion flux of 222 Rn across the sediment-water interface could be calculated from the following equation [69]: where F sed is the diffusion flux of 222 Rn from the bottom sediments (Bq m −2 h −1 ); λ is the decay constant of 222 Rn (7.56 × 10 −3 h −1 ); φ is the porosity of sediments (0.49 for the surface sediments in our study area); C eq and C 0 are the activities of 222 Rn in the pore water and the overlying seawater. Here, we took an average value of 533 Bq m −3 as the activity in pore water. D m is the molecular diffusivity coefficient, which fell into a range between 1.22 × 10 −5 and 1.64 × 10 −5 cm s −1 in our study, through Equation (14). Based on Equation (16) The former is for the flooding tide period; while the latter is for the ebbing tide period. C w is the observed activity of 222 Rn in the coastal water, and C S represents the 222 Rn activity in the offshore seawater (0.28 Bq m −3 ). The term, b, is the return flow factor (i.e., the percentage of the tidal prism that returns from the offshore region during a rising tide). To calculate the return flow factor (b) in the Gaomei Wetland, we employed a tidal prism model, which is based on a mass balance of salinity [70], and obtained b = 0.39, which means offshore seawater occupies more than 60% of the total volume discharging into our study area during flood tide. We calculated the flushing time (T f ) of 0.64 days of the Gaomei Wetland using the following equation [71,72], as well: where V is the total volume of the Gaomei Wetland (1.05 × 10 6 m 3 ); T is the tidal period (12.5 h); P = 1.40 × 10 6 m 3 is the tidal prism, evaluated from P = QT [71,72]. Q is the average water discharge of the Da-Chia River, which is displayed in the next section. The Da-Chia River is one of the 222 Rn sources that we have to consider for obtaining the 222 Rn attributed to SGD. The averaged water discharge (Q River ) of the Da-Chia River was 31 m 3 s −1 [39] and the 222 Rn activity in river water was on average 1186 ± 216 Bq m −3 . The resulting 222 Rn-river flux from the Da-Chia River into the Gaomei Wetland was 18.9 Bq m −2 h −1 .

SGD Rate Estimating
With all known sources and sinks of 222 Rn for the Gaomei Wetland, we estimated the 222 Rn flux attributed to SGD (F SGD ) and then converted this flux to an SGD flux (Q SGD ) by dividing by the 222 Rn end-member concentration in groundwater (C GW : 9187 ± 595 Bq m −3 ). The estimated SGD fluxes into the Gaomei Wetland over the experimental periods during both dry and wet seasons fluctuated with the tidal cycles, showing higher fluxes during the low-tide period, and lower ones during the high-tide period (Figure 8). This indicated that SGD fluxes across the seabed was not in a steady state over the time scales of tidal period, but response to the tides that reached a peak about every 12 h. For the dry season, the SGD fluxes varied between 0.2 and 25 cm d −1 , with an average value of 12.9 ± 7.7 cm d −1 . For the wet season, it ranged from 0.1 to 47 cm d −1 , with a mean value of 19.0 ± 13.6 cm d −1 (Figure 8, Table 3). The slightly high SGD fluxes occurring at spring tide of wet season might have resulted from a much stronger tidal pumping, coupled with a large hydraulic gradient in land. It implied that the general fluctuation of SGD into the Gaomei Wetland responds to both hydraulic gradient and tidal effect. A dry season (i.e., less groundwater recharge) caused a relatively low difference of hydraulic heads between groundwater and seawater, which might result in decreasing SGD fluxes, or even seawater intruding into the coastal aquifer [73]. In the wet season, high groundwater recharge, which implies an increasing hydraulic gradient in land, conversely leading to a rising amount of freshwater from the seawards land flows; thus, fresh groundwater discharging offshore is enhanced [73]. On the other hand, due to tidal pumping, the recirculated seawater mixed with fresh groundwater, moves out through the shallow aquifer and sediments, at low tide. The difference of hydraulic heads between groundwater and seawater increases, and then causes larger SGD fluxes [11]. At high tide, the difference of hydraulic heads between groundwater and seawater decreases, so the SGD fluxes decreased as well, or even seawater flows into the shallow aquifer and sediments [11].

Uncertainty Analysis
It should be pointed out that the estimates of the SGD fluxes through the mass balance model of 222 Rn are always subject to inherent uncertainties associated with its sources and sinks. In our calculation, the uncertainties attributed to the diffusion from sediments, the atmospheric loss, and the fluvial input are less important. A 10% change on 222 Rn fluxes of diffusion from sediments, atmospheric loss, and fluvial input result in 0.5 %, 0.1%, and 5% changes of the SGD flux, respectively. Uncertainty analysis exhibited that the SGD estimates are highly sensitive to the 222 Rn end-member value of groundwater, with a 10% shift in the 222 Rn concentration of groundwater results in 10% change in SGD flux.
To confirm if our estimates are credible, we roughly calculated the SGD rates through Darcy's Law: where V is the flow rate of groundwater (m s −1 ), i is the hydraulic gradient, and k is the hydraulic conductivity (m s −1 ). For groundwater in the catchment of the Da-Chia River, i ranges from 0.007 to 0.01; while k varies between 2.0 × 10 −5 and 1.0 × 10 −3 m s −1 [41]. Thus, the groundwater-flowing rate falls into a range between 1.20 and 86.4 cm d −1 (Table 3). In addition, based on a seepage-meter measurement, Lin [24] reported an SGD flux into the Gaomei Wetland of 1.68 cm d −1 , which consisted of 65% of fresh groundwater. All of these results show good agreements with our SGD estimates from a 222 Rn mass balance model.

SGD-born DIC and Nutrient Fluxes
We used the average SGD flux estimates (12.9 ± 7.7 cm d −1 for dry season; 19.0 ± 13.6 cm d −1 for wet season) and the average DIC concentration in groundwater (1.63 ± 0.09 mM; [59]) to calculate the SGD-born DIC flux to the Gaomei Wetland. The results were 1.5 × 10 6 and 2.2 × 10 6 mol d −1 in the dry and wet seasons, respectively. These fluxes were~26% of the DIC flux from the Da-Chia River (7.1 × 10 6 mole d −1 ). Our SGD-born DIC flux was comparable to the result obtained in the Okatee Estuary, South Carolina, (1.5 × 10 6 mole d −1 ; [72]), where is a similarly scaled estuary, as compared to our study area. Our results were much lower than those obtained in the North South China Sea (NSCS) shelf, off the Perl River (419 × 10 6~9 50 × 10 6 mole d −1 ; [14]). Here, a high DIC concentration in the groundwater led to a high SGD-born DIC flux.
Once exported to the coastal zone, the SGD-induced DIC could be consumed by biological production, which was stimulated by the nutrients supplied from underground water [14]. Lin [24] reported that the SGD-derived nutrient fluxes into the Gaomei Wetland were 53.5 × 10 3 and 1.71 × 10 3 mole d −1 for N and P, respectively. They were about 57% for N and 70% for P, as compared to those exported from the Da-Chia River. In addition, the N to P ratio was 31.3, well above the Redfield ratio of 16 [74], which suggests that SGD-derived nutrient inputs might drive the Gaomei Wetland to P limitation. Based on the Redfield ratio (C:N:P = 106:16:1; [74]), the SGD-derived nutrient might be responsible for a new production up to 0.18 × 10 6 mole C d −1 , which meant that~11% of SGD-derived DIC would contribute to the new production. This was comparable to the result observed in NSCS shelf, off the estuary of the Pearl River [14]. We confirmed that SGD is an important pathway for carrying nutrients and DIC from land to the coastal ocean and plays a non-negligible role in the biogeochemical cycle.

Summary
We reached the following summaries based on a time-series investigation of 222 Rn, 224 Ra ex , 228 Ra, δ 18 O and δD in coastal water, and pore water over 2 tidal cycles in both dry (9~10 May 2014) and wet (27~28 August 2014; encountering a spring tide) seasons, at a station 500 m offshore in the Gaomei Wetland, Western Taiwan.
A. In the wet season, 224 Ra ex and 228 Ra activities in the coastal waters were high (>20 dpm 100 L −1 ) during the low-tide period, but relatively low (<10 dpm 100 L −1 ) during the high-tide one; whereas an opposite trend was displayed in the dry season. The much higher 224 Ra ex and 228 Ra activities in the low-tide water in the wet season might not have resulted from the contribution of Ra desorption from marine suspended particles; coupled with the δ 18 O and δD signatures of water, it indicated a source of RSGD. B. The time-series observations of 222 Rn showed the same trends regardless of wet and dry seasons, with higher activities at low tides and lower ones at high tides. Based on a 222 Rn mass balance model taking all sources and sinks into account, we estimated an SGD flux into the Gaomei Wetland, ranging from 0.1 to 47 cm d −1 , which varied with tide fluctuation, with low rate during the high-tide period, but a high one during the low-tide period. This confirmed the importance of tidal pumping in driving SGD. The slightly higher fluxes during the wet season (average SGD flux = 19.0 ± 13.6 cm d −1 ) compared to the dry season (average SGD flux = 12.9 ± 7.7 cm d −1 ) indicated a stronger tidal pumping coupled with a large hydraulic gradient in land. C. The DIC fluxes via SGD into the Gaomei Wetland were estimated to be 1.5 × 10 6 and 2.2 × 10 6 mol d −1 in the dry and wet seasons, respectively,~26% of the river-borne DIC flux from the Da-Chia River. In addition, the SGD-derived nutrient fluxes were 53.5 × 10 3 and 1.71 × 10 3 mole d −1 for N and P, respectively, which were about 57% for N and 70% for P, as compared to those exported from the Da-Chia River. Such exports of nutrients and the DIC fluxes via SGD could have had an impact on the coastal biogeochemistry in the Gaomei Wetland. Funding: This research received no external funding.