Seasonal and Spatial Production Patterns of Dissolved Inorganic Carbon and Total Alkalinity in a Shallow Beach Aquifer

Dissolved inorganic carbon (DIC) and total alkalinity (AlkT) fluxes to the nearshore ocean can directly impact the rates of primary production, coral reef formation, coastal ocean acidification, and continental shelf ecology. Current understanding of the transformations that DIC and AlkT undergo as they move from land to sea are limited, leading to difficulties in estimating future DIC and AlkT export that may be altered under a changing climate. While much research has focused on carbon fluxes in carbon-rich mangroves and coastal wetlands, DIC and AlkT transformations and distributions in sandy beach aquifers, which are comparatively carbon-poor, have not been studied as extensively. We monitored DIC and AlkT concentrations in a sandy beach system over six sampling events spanning two years. Substantial changes to DIC and AlkT occurred along subsurface flowpaths due to aerobic respiration and anoxic reactions, resulting in an additional mean flux to the ocean of 191 and 134 mmol/d per meter length of shoreline, respectively. The chemical alterations occurred within the saltwater-freshwater mixing zone beneath the beach surface. Both aerobic and anaerobic reactions actively contributed to DIC and AlkT production within the system, as indicated by DIC: AlkT and dDIC:dAlkT ratios relative to the theoretical dilution line. The work indicates that beach aquifers support active transformation of inorganic carbon and highlights a potentially important and overlooked source of DIC and AlkT to coastal systems.


INTRODUCTION
Controls of dissolved inorganic carbon (DIC) and total alkalinity (Alk T ) fluxes to estuaries and the ocean are gaining attention as increased carbon emissions, changing land use, and climatic feedbacks alter the amount of CO 2 entering the atmosphere, contributing to climate change (Bauer et al., 2013;Regnier et al., 2013). There is debate concerning the role of coastal zones in modulating the inorganic carbon cycle within the context of climate change. While continental shelves take up atmospheric CO 2 via primary production, estuaries, especially in lower latitudes, have been shown to deliver large amounts of CO 2 to the atmosphere via inorganic carbon fluxes from groundwater, wetlands, and rivers (Jiang et al., 2008). The relative fluxes of reactive organic carbon (OC), DIC, and Alk T within and across these systems is likely to play a fundamental role in the rates of ocean acidification and climate evolution (Cai, 2011;Hu and Cai, 2011;Liu et al., 2021). Coastal ecosystems are key components of the global carbon cycle owing to the disproportionately large amount of biogeochemical activity, including carbon cycling, relative to surface area (Cai, 2011). While the magnitude of organic and inorganic carbon fluxes from mangroves, estuaries, volcanic aquifers, and coral reef systems have been recently studied (Cyronak et al., 2013;Cardenas et al., 2020;Connolly et al., 2020), the role of beach aquifers in carbon fluxes to the coastal zone has not been thoroughly assessed. Within beach sediments, the magnitude of production and consumption of DIC and Alk T have been shown to be site-specific and vary seasonally (Chaillou et al., 2014;Liu et al., 2021). These spatial and temporal variations in DIC and Alk T production and consumption are due to the complex interactions between freshwater and seawater, which can occur at various scales within beach systems (Michael et al., 2005), and the lack of labile organic carbon compared to carbon-rich systems such as mangroves or wetlands.
Groundwater discharging across the seabed (submarine groundwater discharge; SGD) consists of fresh, brackish, and saline groundwater (Robinson et al., 1998;Taniguchi et al., 2019;Santos et al., 2021). Nutrient concentrations in SGD are often higher than concentrations in riverine surface water, thus SGD fluxes can negatively affect coastal water quality, contribute to eutrophication, and compromise coastal marine ecosystems (Johannes, 1980;Moore, 1999;Slomp and Van Cappellen, 2004;Guo et al., 2020). Further, SGD and associated chemical fluxes are biogeochemically dynamic at tidal to seasonal time scales (Beck et al., 2008;Santos et al., 2009;Kim et al., 2017;Liu et al., 2018;Waska et al., 2019). Prior to discharge, dissolved chemicals such as reduced iron, nitrate, and ammonium are transported through the beach aquifer, where active carbon cycling stimulates redox processes that modulate chemical fluxes from land to the sea (Cai et al., 2003;Beck et al., 2017;Kim et al., 2017;Kim et al., 2019). Recent studies have demonstrated that these redox processes affect the net production and sequestration of CO 2 in intertidal porewater (Liu et al., 2021). Globally, coastal groundwater fluxes of HCO − 3 to the coastal zone are estimated to be 32-351% of fluxes derived from weathering of silicate rocks (Pain et al., 2019). The uncertainty in this estimate reflects the need to more accurately quantify groundwater DIC and Alk T fluxes to the ocean.
In some sandy beaches, a persistent and dynamic freshwaterseawater (FW-SW) mixing zone forms in the intertidal zone (Michael et al., 2005;Greskowiak, 2014). Wave and tidal action transport seawater across the sandy beachface. The seawater infiltrates the porous beach sediments and mixes with seawardflowing fresh groundwater, creating a wave-and tidally-driven cell of circulating seawater. Because seawater contains high amounts of labile organic carbon, infiltrating seawater can support biogeochemical transformations in the beach, such as aerobic respiration, denitrification, iron oxidation/reduction, and sulfate reduction (Charette and Sholkovitz, 2002;McAllister et al., 2015;Reckhardt et al., 2015;Beck et al., 2017). Thus, the circulation zone hosts reactions that alter chemical concentrations and fluxes of various reactive dissolved species, including oxygen, organic carbon, and nitrate (Ullman et al., 2003;Hays and Ullman, 2007;McAllister et al., 2015;Reckhardt et al., 2015;Couturier et al., 2017;Kim et al., 2017;Y. Liu et al., 2017). Geologic heterogeneity within beach and tidal riparian aquifers can influence the spatial distribution of the reactions as well as the chemical fluxes in discharging groundwater Wallace et al., 2020).
Few studies have observed changes in DIC and Alk T concentrations and fluxes to the ocean due to biogeochemical reactions in shallow nearshore sediments (Burt, 1993;Cai and Wang, 1998;Lee et al., 2006;Cai et al., 2015). Cai et al. (2003) showed that dissolved species concentrations (DIC, Alk T , pH, pCO 2 ) in the saline-fresh groundwater mixing zone at the North Inlet, South Carolina, USA were altered by carbonate dissolution and organic carbon degradation via sulfate-reduction. Chaillou et al. (2014) similarly estimated changes to DIC and Alk T flux across a shallow SGD zone in the Gulf of St. Lawrence, Quebec, Canada. The estimates were obtained by multiplying the fresh groundwater flux by the concentration of DIC and Alk T in shallow beach groundwater (10 cm beneath the beach surface). Kim et al. (2017) showed that aerobic respiration and CO 2 production along the freshwater-saline water mixing zone enriched pCO 2 and DIC concentrations in the mixed region. These studies suggest that sandy beach aquifers can serve as local sources, and sometimes as sinks, of DIC and Alk T to the adjacent coastal system. Recently, Liu et al. (2021) documented the temporal variations of DIC and Alk T distributions and production within an intertidal aquifer in Hong Kong, and suggested that the capacity of beach to serve as a DIC and Alk T source or sink is site-specific.
The goal of this study was to explore the seasonal variations to DIC and Alk T production due to aerobic and anaerobic respiration in a sandy beach system. To do this, we quantify DIC and Alk T production over two seasonal timescales to provide insight into how beaches moderate DIC and Alk T export to coastal waters. We present the results from six cross-sectional porewater sampling events spanning two years at a sandy beach. The results show that both aerobic and anaerobic biogeochemical transformations in sandy beaches can an important control on groundwater DIC and Alk T export to the coastal ocean, with implications for regional to global inorganic carbon budgets.

FIELD SITE AND METHODS
Field Site range = 1.42 m), sandy (median grain size = 540-668 µm) beach with a slope of 0.1 and is protected from storm erosion by two offshore breakwaters. The upland hydraulic head varies seasonally and in response to rainfall, and was monitored in a well 33 m landward of the dune using a water level logger programmed to record at 15-minute intervals. The recorded water levels were used to determine seasonal changes to the freshwater hydraulic gradient over approximately a year and half. Precipitation data were obtained from the National Oceanic and Atmospheric Administration (NOAA) Georgetown Sussex Country Airport Station USW00013764, located 26.8 km from the study site.
Previous research at the Cape Shores site has confirmed the persistence of a perennial, but seasonally mobile, circulation cell within a zone of active mixing between seawater and fresh groundwater (Ullman et al., 2003;Heiss and Michael, 2014;Kim et al., 2017;Kim et al., 2019). The circulation cell at this site hosts transformations of various chemical constituents, such as nitrate, ammonium, iron, and carbon (Ullman et al., 2003;Hays and Ullman, 2007;McAllister et al., 2015;Kim et al., 2017). At this site, aerobic respiration has been found to be the predominant pathway affecting inorganic carbon, with minimal effects from carbonate dissolution (Ullman et al., 2003;Kim et al., 2017).

Methods
Porewater was sampled using multi-level samplers constructed from ½'' PVC pipe with polyethylene tubing extending down to different sampling depths (Heiss and Michael, 2014;Kim et al., 2017;Kim et al., 2019). On average, the samplers were~3 m long with polyethylene sampling ports spaced 75 cm apart. The sampling port tips were covered with mesh to prevent clogging. The samplers were installed along a shoreperpendicular transect across the beachface, extending 167 meters seaward of the dune, across the beachface mixing zone, and to a depth of~4 meters below the beach surface. Sampling of all parameters and porewater collection and preservation took 6-8 hours, and measurements therefore represent different times in a tidal cycle. However, previous research at this site has shown that intertidal salinity variability over semidiurnal cycles is minimal (Heiss and Michael, 2014); the spatial patterns presented here can be taken to be diurnally-averaged snapshots.
A peristaltic pump with an attached filtration cartridge (0.7 µm glass fiber filter) was used to collect water samples. Salinity, dissolved oxygen (DO) saturation (%), and pH were measured using a YSI Professional Plus meter (YSI Incorporated; Yellow Spring, Ohio, USA) with an attached flow-through chamber. Porewater was pumped directly from the sampling ports through the filter and into clean and borosilicate bottles pre-rinsed with filtered water. The ports had a vertical separation distance of 50-80 cm and the volume of water pumped on average for sampling (1-3 liters) was limited to prevent vertical exchange of porewater across sample ports. Samples were immediately fixed with 80 µl saturated HgCl 2 , and stored in the dark at 4°C until analysis. DIC was measured in the laboratory as CO 2 following acidification of using an Apollo SciTech AS-C3 DIC Analyzer. Alk T was determined by Gran titration (Gran, 1952) using the Apollo SciTech AS-ALK2 Alkalinity Titrator (see Huang et al., 2012;Kim et al., 2017;Kim and Heiss, 2021 for more details on methodology). The analytical precision for three replicate samples was ± 2-3 µmol kg -1 for both DIC and Alk T. To obtain seasonal cross-sectional snapshots of the measured parameters, samples were collected at each sample port six times over a two-year period ( Figure 1).
DIC and Alk T production at each sample port was calculated as the difference between the measured DIC or Alk T concentration and the concentration at the theoretical dilution line (see Methods for details). The average produced concentration was multiplied by the brackish SGD rate to estimate excess DIC and Alk T fluxes to the coastal ocean at each sampling event. Brackish SGD for each sampling event was obtained from a numerical model of the Cape Shores beach aquifer (calibrated to the period May 2012 -May 2013; Heiss and Michael, 2014). The brackish SGD rate on the day and time of each sampling event was extracted from the continuous 1-year model by summing SGD rates at model cells along the aquiferocean interface with a salinity within the range measured in the field during that month. The same calculation was performed for all months with the exception of November 2014 where salinity in the circulation cell was very high, suggesting the production zone extended outside of the measurement area. To compensate for production outside of the measurement area for November 2014, the salinity range across all months was used to define the brackish flux component for that month. For all brackish fluxes, density-driven brackish SGD from the deeper saltwaterfreshwater interface was excluded so that only brackish fluxes from the intertidal circulation cell where measurements were collected were included.

Spatial Distribution of Salinity, Oxygen, DIC, Alk T , and pH
The salinity distribution in the intertidal zone varied with seasonal changes to the freshwater hydraulic gradient and tidal amplitude. Groundwater was freshest near the most landward sampler and progressively increased in salinity as it mixed with infiltrating seawater within the intertidal zone ( Figure 2; salinity panels, zone of seawater infiltration). This created a persistent, but positionally variable, "intertidal circulation cell", defined by salinity. Salinity was generally highest below the high tide line where the seawater infiltration rate is highest (Heiss and Michael, 2014;Heiss et al., 2015). Fresh to brackish groundwater discharged to the ocean seaward of the center of the intertidal circulation cell, consistent with modeled salinity distributions at this site and in other beach aquifers (Robinson et al., 2007;Abarca et al., 2013;Beck et al., 2017;Kim et al., 2017). As freshwater flux decreased heading into the winter months, the fresh-brackish discharge zone was displaced seaward as the intertidal circulation cell expanded due to greater seawater infiltration across the upper beachface (e.g., June 2014 to November 2014 and May 2015 to September 2015). The extent of the DO saturation followed the lateral extent of salinity, but was quickly depleted at depth once infiltrated, leading to lower saturation values than would be expected with conservative mixing (Figure 2, DO panels). As wave runup is saturated in oxygen due to its contact with the atmosphere, this indicates active oxygen consumption via aerobic respiration within the beach sediments. Rapid oxygen consumption creates anoxic zones suitable for anaerobic reactions in the deeper regions of the intertidal zone (Charbonnier et al., 2013;Kim et al., 2017).
The spatial distributions of DIC and Alk T generally followed salinity, with higher DIC and Alk T concentrations correlating with higher salinities. This led to a wide range of DIC and Alk T concentrations within the intertidal zone over the six sampled months, with DIC spanning 423-2378 µmol/kg and Alk T ranging from 423-2405 µmol/kg (Figure 2). DIC and Alk T concentrations also displayed seasonal patterns, similar to the temporal variability in salinity through sample period spanning two years. In the freshest regions of the aquifer across all 6 sampling events, DIC ranged from 423-1427 µmol/kg and Alk T ranged from 423-1108 µmol/kg. Interestingly, the highest DIC and Alk T concentrations did not always correspond to the locations with highest salinity (Figure 2, stars), indicating nonconservative behavior and active production of excess DIC and Alk T within the intertidal zone. Further, while Alk T is greater than DIC in seawater (Cai et al., 2003), Figure 2 shows that DIC concentrations in the porewater were higher than Alk T (note that DIC and Alk T colorbar scales are the same and DIC panel colors are more saturated). This is due to the production of DIC within the beach sands prior to its coastal discharge. In some samples, Alk T was artificially low because of HgS precipitation during HgCl 2 fixation; thus, pH could not be calculated using the CO2SYS tool as a balance between DIC and Alk T . Instead, results from the YSI sensor are presented. pH is high (> 8.0) in the freshwater dominated regions (except September 2014) but reduces quickly to less than 7.0 in the freshwater-salt mixing zone, where DIC is produced (see Methods). Such deviations in pH were observed in all sampling events, suggesting active transformation of inorganic carbon within the beach system.

Production of DIC and Alk T
DIC and Alk T production (i.e. excess DIC and Alk T , only DIC presented) was calculated for each porewater sample by subtracting the measured concentrations from the expected concentrations assuming conservative mixing of the freshest and most saline samples collected on that date. DIC and salinity was used to select the fresh and saline end-members, as DIC alterations from aerobic alterations were expected prior to any Alk T alterations. For each month, a linear slope from the highest measured DIC and Alk T concentration to the fresh endmember was calculated as the theoretical dilution line (TDL; Figure 3). We used the sample with the lowest salinity as a proxy for the freshwater endmember. The sample with the lowest salinity was often among the lowest-DIC concentration samples, supporting its selection as the fresh endmember prior to DIC-altering reactions. One exception was July 2014, where the freshwater endmember sample at salinity 0.1 had 300 µmol/kg more DIC than a neighboring sample (1.5 m deeper) that had a salinity of 6. This is potentially due to thermal effects, as colder temperatures in the freshwater endmember would increase DIC concentrations and make samples fall below the TDL. Nonetheless, this approach would generate a conservative estimation of internal DIC production, as the true endmember is likely lower, for example, in the case of November 2014 ( Figure 3). The DIC concentrations of the fresh end-members were also similar to the Delaware River fresh end-member samples at the Trenton river, located up-stream from the Cape Shores site (Joesoef et al., 2017& unpublished data from A. Joesoef, 2017, which further supports FIGURE 3 | Salinity vs. DIC. The theoretical dilution line assuming freshwater DIC = 0 is shown as a solid line, and assuming a freshwater DIC endmember from the lowestsalinity sample shown as a dotted line. Greatest production was measured in samples within the discharging water (pink symbols). See Supporting Information for sample grouping method that utilized cell geometry and salinity (November 2014 prioritized cell geometry in sample grouping due to a broad saline region). the chosen endmembers. Average salinity of the bay water was 28. While seawater DIC and Alk T for seawater was not measured for the sampling dates, samples within the infiltration zone ( Figure 3; green triangles) were close to seawater salinity and taken as saline endmembers prior to chemical transformation.
Production of DIC and Alk T was spatially correlated to the freshwater-saltwater mixing zone (Figure 4). The greatest deviations (and therefore, the greatest DIC and Alk T production) were observed in samples that were within the brackish mixed zone or the discharge zone within the circulation cell (See Supplementary Material). The magnitude of highest production across all months varied from 762-1746 µmol/kg for DIC and from 535-1523 µmol/kg for Alk T , indicating substantial year-round production of DIC and Alk T within the intertidal aquifer (Tables 1, 2). For flux calculations, the average deviation across all samples were taken to and multiplied with the brackish SGD flux. This active production led to on average an additional 191 mmol/d DIC and 134 mmol/d Alk T flux to the ocean per meter of shoreline (Tables 1, 2).
The only exception was July 2014, when the freshwater endmember sample selected had slightly elevated DIC concentrations compared to other neighboring samples (salinity 2-6). Nonetheless, the DIC concentrations for discharging water samples (Figure 3, pink symbols) were higher than TDL, indicating a presence of DIC-altering reactions such as aerobic respiration (also shown in Figure 4 in Discussion). Thermal effects seem probable for the slightly elevated DIC for the freshwater endmember as carbonate dissolution has not been observed at this field site as a mechanism for DIC removal, where quartz is the dominant composition of sand (Kim et al., 2019). The negative flux of Alk T in July 2014 can be attributed to the active FeS precipitation that has been observed at this site (McAllister et al., 2015), that may have been more active during this month compare to other months. Further discussion on the diagenetic pathways are discussed in Section 3. Nonetheless, the seasonality observed in the altered fluxes for DIC and Alk T are in line with the behavior observed by Liu et al. (2021) in Tolo Harbor, Hong Kong.
The changes to DIC and Alk T flux occurred regardless of season, and was unrelated to seasonal average porewater temperature or inland freshwater head ( Table 1). These changes DIC and Alk T fluxes may alter the carbon chemistry of coastal waters through a decrease of pH, causing acidification.

Spatial Characteristics of DIC and Alk T Production
Spatially, production of DIC and Alk T was greatest along the mixing zones between freshwater and saline water within the circulation cell. Table 1 shows the salinity of samples that were chosen to calculate maximum excess DIC and Alk T production for each month. The highest concentrations of excess DIC and Alk T most often occurred in brackish mixing zones, save for November 2014 that consistently had higher salinities due to an overall more saline system. For all months, the highest excess DIC/Alk T was observed in brackish or near-discharge zones.
The cross-sections of excess DIC and Alk T (Figure 4) demonstrate the heterogeneous distribution of the reactions in driving DIC and Alk T production. DIC and Alk T production were the highest in the mixing zone between freshwater and saline seawater along the perimeter of the circulation cell and was spatially discontinuous, occurring below both the upper and FIGURE 4 | Spatial distribution of excess DIC and T-Alk, and DIC : Alk T ratios over seasonal timescales. Salinity 15, 20, and 25 contours are shown. Higher ratios are often found along freshwater-seawater mixing areas due to aerobic respiration, and subsequent anaerobic respiration supported by oxidants in freshwater. High DIC : Alk T ratios near the discharge zone (seaward samplers near the 160 m mark) are due to high sulfate reduction activity. Solid black line indicates topography, with mean high-tide and mean sea-level indicated with dotted horizontal lines. lower beachface (Figure 4). This is consistent with previous studies that demonstrated that aerobic respiration of marinederived organic matter in intertidal sediments is highest along subsurface mixing zones and that CO 2 accumulates along flowpaths leading towards the discharge zone (Charbonnier et al., 2013;Seidel et al., 2015;Kim et al., 2017;Kim et al., 2019). The high production rates in the mixing zone beneath the upper beachface can be attributed to aerobic respiration supported by an abundant supply DOC and DO in infiltrating seawater, while the activity beneath the lower beachface is likely due to anaerobic respiration, which subsequently occurs along the circulating and progressively anoxic flowpaths (Figure 4, pink arrow in July 2014 Excess DIC panel). The anaerobic respiration may be supported by oxidants in upwelling freshwater (e.g., nitrate).

Diagenetic Pathways of DIC and Alk T Production at Cape Shores
The diagenetic process that alters DIC and Alk T within the beach system can be revealed by DIC : Alk T ratios and the slope of dDIC:dAlk T . Without major biogeochemical transformations, we expect the DIC : Alk T ratio to fall on 1 (Cai et al., 2003). The observed DIC : Alk T ratios for both fresh and saline endmember samples for most months fell on close to the 1:1 line ( Figure 5, blue and orange dots, respectively), indicating that a 1:1 ratio can be expected for samples prior to undergoing biogeochemical reactions that would alter DIC or Alk T concentrations and ratios. DIC : Alk T ratios near 1:1 in beach aquifer systems and nearshore seawater have similarly been observed at other sites (e.g. Chaillou et al., 2014;Liu et al., 2021). The low oxygen saturation within the system indicates elevated levels of aerobic respiration (Figure 2). Because aerobic respiration only influences DIC production through the production of CO 2 (Cai et al., 2003), it is expected that DIC : Alk T/ dDIC:dAlk T ratios will reach near infinity. However, while there is a slight positive deviation in DIC : Alk T and dDIC: dAlk T ratios for non-endmember samples, most of the samples fall along the 1:1 line (Figures 5, 6). This indicates a strong presence of anerobic reactions in addition to the aerobic reactivity. Under anoxic conditions, alterations to alkalinity  The salinity column indicates the salinity of the highest DIC sample during that month, which was never the highest-salinity sample. Average produced DIC was multiplied with the brackish (salinity of the highest DIC sample) SGD flux to yield changes to DIC flux via reactions into the oceanic system. For simplicity, SGD was assumed to be 1013kg/m 3 , the mean value between freshwater and saline bay water.
can occur via denitrification, sulfate reduction, and sulfur-iron interactions. Previous research at this study site has also shown that organic matter degradation below the lower beachface proceeds via denitrification and sulfate reduction rather than aerobic respiration (McAllister et al., 2015;Kim et al., 2017;Kim et al., 2019). Carbonate dissolution rates at Cape Shores are expected to be minimal due to its quartz-dominated sedimentology. Denitrification and sulfate reduction, however, are expected to increase Alk T and DIC in a ratio close to 1 via H + consumption (Canfield et al., 1993;Cai et al., 2003;Cai et al., 2010a;Cai et al., 2010b): H + + 2CH 2 O+SO 2− 4 ! 2CO 2 + HS − + 2H 2 O While FeO(OH) 2 reduction could serve as a source of Alk (Van Cappellen and Wang 1996;Cai et al., 2010a;Cai et al., 2010b), we expect that such Alk production can largely be cancelled out because precipitation of FeS from the diagenetic production of HSin the presence of Fe 2+ also consumes two moles of alkalinity per mole of HSremoved (Canfield et al., 1993;Cai et al., 2003;Cai et al., 2010a;Cai et al., 2010b;McAllister et al., 2015): Fe 2+ + HS − ! FeS + H + Sulfate reduction, however, has been observed at Cape Shores and maintains DIC : Alk T and dDIC:dAlk T ratios near 1 (Cai et al., 2010a;Cai et al., 2010b;Su et al., 2021). Although Fe 2+ , SO 2− 4 and HSwas measured for July 2015 only (Supporting Information, Figure S1), high Fe 2+ concentrations and a sulfate deficit (despite the high salinity) in the discharge zone is supported by previous research (McAllister et al., 2015). Another indicator of high S presence was the precipitation of HgS (thereby decreasing measured Alk T ; Figure 5, red squares) during HgCl 2 fixation when the samples were being preserved.
The reaction stoichiometry discussed above indicates that the Cape Shores system has both active aerobic and anaerobic reactions that produce excess DIC and Alk T at rates that preserve the ratio between them near 1. The biogeochemical processes that alter DIC and Alk T concentrations along circulating flow paths (indicated by salinity; an example flowpath is given by a pink arrow in Figure 4, July 2014 Excess DIC panel) can be summarized as follows: Infiltration of DOC-and DO-rich seawater across the upper beachface supports high rates of aerobic respiration below the high tide line, increasing DIC relative to Alk T . Farther along the circulating flow path, limited DO availability leads to anoxic conditions and conditions favorable to anoxic reactions, including denitrification, supported by nitrate transported into the beach by inflowing freshwater, and sulfate reduction towards the discharge zone. Both denitrification and sulfate reduction produce DIC and Alk T at equal rates, resulting in DIC : Alk T ratios near 1. While we did not evaluate the contribution of aerobic and anaerobic reactions to DIC and Alk T chemistry, the findings demonstrate that freshwater-saltwater mixing zones within intertidal beaches strongly affects the inorganic carbon chemistry of groundwater discharging to the ocean, increasing concentrations up to fourfold.

Potential Controls and Effects of DIC and Alk T Production at Cape Shores
Excess DIC and Alk T production and export to surface water is likely influenced by a number of factors that affect the biogeochemical reactivity of beaches. In this work, DIC and Alk T, while persistent seasonally, were not strongly linked to seasonal variations in the inland water table or average porewater temperature ( Table 1). The average flux of 191 mmol-m -1 -d -1 DIC across seasons was lower than observations in other beach systems, such as the 2-3 mol-m -1 -d -1 DIC flux observed at the Magdalen Islands in the Gulf of St. Lawrence (Chaillou et al., 2014) and the 51-93 mol-m -1 -d -1 DIC flux observed at the Tolo Harbor system (Liu et al., 2021). The 134 mmol-m -1 -d -1 Alk T for Cape Shores, however, was higher than the Magdalen Islands system (2-8 mmol-m -1 -d -1 Alk flux) and lower than the Tolo Harbor system (48-73 mol-m -1 -d -1 Alk flux). This highlights the local difference in biogeochemical reactivity of beaches and may be controlled by the size and extent of the intertidal zone.
Generally, DIC and Alk T production will depend on oxygen and other oxidant (e.g., nitrate, sulfate, iron) availability, organic carbon availability and reactivity, in addition to porewater temperature. Seasonal variability of ocean primary production would likely control the amount of carbon introduced into the beach subsurface, and thus the type and extent of reactions capable of modifying DIC and Alk T concentrations. Conditions that lead to rapid oxygen consumption, such as warm coastal porewater temperatures (Cogswell and Heiss, 2021) or high labile dissolved organic carbon availability (Spiteri et al., 2008), would also control DIC and Alk T production ratios through their control on the redox environment. Additionally, influx zones of labile particulate organic carbon in infiltrating seawater can be offset relative to carbon-fueled reactions within the aquifer due to retardation during transport of particulate organic carbon, as well as leaching (Kim et al., 2020). DIC and Alk T are likely to be affected by overall beach carbon availability, including less mobile particulate organic carbon pools, and therefore may explain the lack of seasonal trends in production. Nonetheless, this work highlights seasonally active alterations to DIC and Alk T and shows that despite low in-situ organic carbon availability relative to other systems such as mangroves or wetlands, sandy beaches can contribute to nearshore pH and carbon chemistry. In more tropical regions, this added flux of DIC and Alk T could impact coral ecosystems via seawater acidification or coral calcification. Such implications would require site-specific analysis and in-depth monitoring of DIC and Alk T concentrations and fluxes within tropical coastal beaches.

CONCLUSIONS
This study demonstrates that mixing between freshwater and seawater in beach aquifers can support substantial production of DIC and Alk T through redox processes. Porewater salinity in the intertidal zone varied over a two-year sample period, driven by changes in terrestrial hydraulic head and tidal conditions. Dissolved oxygen was high in the shallow region of the aquifer and was depleted at depth, indicating active aerobic respiration and oxygen consumption. Aerobic respiration led to DIC production hotspots in the freshwater-seawater mixing zone, resulting in DIC : Alk T ratios slightly higher than 1. While in the absence of anaerobic reactions DIC : Alk T ratios would reach infinity, excess Alk T within the mixing zone also indicated anoxic reactions, such as denitrification and sulfate reduction, which maintained near 1:1 DIC : Alk T ratios. The intertidal subsurface at Cape Shores produced excess DIC and Alk T in quantities that varied across the two-year sampling campaign, with a near fourfold increase in DIC and Alk T fluxes to the bay. The results have implications for local to global carbon cycling and carbonate budgets in nearshore regions and highlight beaches as a potentially overlooked system affecting DIC and Alk T fluxes to coastal ecosystems.

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

AUTHOR CONTRIBUTIONS
KK conducted the fieldwork, processed samples, did data analysis, drafted the initial manuscript, and created the figures. JH assisted with fieldwork, assisted with project formulation and provided guidance to data analysis and manuscript writing. HM provided the funds for the research, guided the scientific idea, and edited the manuscript. WU also provided the funds for the research, guided the scientific idea, and edited the manuscript. W-JC provided lab space and instruments, guided the scientific idea, and edited the manuscript. All authors contributed to the article and approved the submitted version.

FUNDING
This research was funded by NSF EAR-1246554 (to HM and WU).