Tidal Marsh Record of Nutrient Loadings in Barnegat Bay, New Jersey

Nutrients and organic compounds in sediments are derived from natural and anthropogenic sources. These sources in estuaries include sediment carried by weathered earth material, atmospheric deposition, tidal currents, agricultural and urban runoff, and sewage discharge. Watershed transports of nutrients and carbon loadings to estuaries vary with inflow volume, whereas oceanic transports vary with frequency of tidal flooding. Distinction between C3/C4 plant inputs identified by δ¹³C and δ¹⁵N signatures in marine phytoplankton help to understand terrestrial or marine source in estuaries (Bianchi, 2007; Peterson and Fry, 1989). Sediments deposited in estuarine tidal marshes provide an excellent means for documenting long-term (e.g., at decadal scales) changes in land use, nutrient, and contaminant loadings. The particle-bound nitrogen (N), phosphorus (P), and carbon (C) compounds in tidal marshes have provided burial history under different environmental conditions (Church et al., 2006; Velinsky, Charles, and Ashley, 2007). Although sediments undergo postdepositional diagenetic remobilization, chemical transformations, biological mixing, and hydraulic processes, the sediment column can record a chronology of accretion rate and nutrient burial in estuaries.

Coring information can be used for the construction of sediment budgets (Schubel and Hirschberg, 1977) and to understand nutrient burial in aquatic environments (Church et al., 2006; Velinsky, Charles, and Ashley, 2007; Velinsky, Sommerfield, and Charles, 2010a,b; Velinsky et al., 2011). Dated tidal marsh sediment cores can provide information on chemical loadings over time (Santschi et al., 2001) and provide an estimate of ecosystem services due to sediment burial. However, specific chemicals, such as nitrogen and phosphorus, can undergo substantial diagenetic remobilization depending on the redox conditions (Burdige, 2006), and as such, accumulation records derived from sediment cores may be limited in certain environments. Additionally, in most of the environments, fidelity of sedimentary records is limited by stratigraphic incompleteness.

River inflow transports phosphorus to the estuaries and oceans, with 90% transported in particulate forms (Froelich et al., 1982). Upon burial, sediment-bound phosphorus has the potential to be recycled back into the water column, depending on factors such as temperature, pH, and oxygen concentration (Froelich et al., 1982). Sediment accumulation and loadings and remobilization rates control retention of buried phosphorus in the sediments. Past studies using the sediment core have identified P burial over time (Church et al., 2006; Kahn and Brush, 1994).

Nitrogen transports mostly in the dissolved form, and its recycling mechanism is much more complex (Castro et al., 2003; Van Breemen et al., 2002). As such, retention of nitrogen in coastal and riverine systems is potentially more restricted than for phosphorus. The stable isotopes of nitrogen (¹⁴/¹⁵N) can help determine the source, fate, and cycling of nitrogen in a waterbody and could be reflected in the nitrogen burial in sediments (Kendall, 1998). In the Delaware Estuary, Church et al. (2006) identified an increase in sediment δ¹⁵N through time as an indicator of nitrogen input because of urbanization. This increase was similar to increases in other estuarine systems with urbanized watersheds (Kendall, 1998; McClelland, Valiela, and Michener, 1997; Ulseth and Hershey, 2005).

The objective of this study was to establish the chronology of carbon, nitrogen, and phosphorus burial in tidal marshes of Barnegat Bay, New Jersey, for insight on bay nutrient loadings since about 1900. Tidal wetlands can be an important repository of bioactive elements, and it is important to understand sediment burial as a removal term in the nutrient budget of the Bay. To meet this objective, the chemical characteristics of sediment cores from marsh depositional areas fringing the barrier island lagoon system of the Bay were analyzed.

The Barnegat Bay (BB)–Little Egg Harbor estuary is located along the central New Jersey coastline in the Atlantic Coastal Plain province (Figure 1). Barnegat Bay is a back-barrier lagoon-type estuary that extends from Point Pleasant south to Little Egg Inlet. The variety of highly productive shallow water and adjacent upland habitats found in this system include barrier beach and dune, submerged aquatic vegetation beds, intertidal sand and mudflats, salt marsh islands, fringing tidal salt marshes, freshwater tidal marsh, and palustrine swamps.

The Barnegat Bay system is composed of three shallow bays (Barnegat Bay, Manahawkin Bay, and Little Egg Harbor), is approximately 70 km in length, 2–6 km wide, and up to 7 m deep. The Bay watershed has an area of nearly 1700 km² and has been extensively developed since 1940. Tidal waters cover approximately 280 km² with a ratio of watershed area to water area of 6.1. In the past two decades, urban land cover has increased with the loss of forested land (Lathrop and Haag, 2007; Paudel et al., 2016). Population of the watershed has increased substantially from the 1940s (40,000) to more than 570,000 year-round resident identified in a recent census. During the summer season, the population can rise to approximately 1,000,000.

In the BB watershed, urban land area has increased since 1986 by 9% (23% urban land area in 1986 vs. 32% in 2012; BBP, 2016). The increase in urban land area was associated with a decrease in forests, farms, and wetlands. Specifically, between 1995 and 2012, an average of ∼4.3 km⁻² of forests, farms, and wetlands were converted to urban land area annually (BBP, 2016). In the same time period, urban population of Ocean County, New Jersey, increased from approximately 471,243 in 1995 to more than 580,945 in 2012 (U.S. Census Bureau, 2016). These changes have intensified nutrient load transport to the Bay over time (Baker et al., 2014).

The water quality of the Barnegat Bay is affected by agricultural runoff and storm water discharge and is influenced by restricted tidal flushing (BBNEP, 2005; Kennish et al., 2007). Approximately 50%–66% of the nutrient load is from surface waters runoff, a substantial amount from atmospheric deposition (22%–40%), and lesser amounts from groundwater inflow (∼10%) (Baker et al., 2014; Bowen et al., 2007; Hunchak-Kariouk and Nicholson, 2001). Kennish et al. (2007), using the NLOAD model framework, estimated nitrogen input of 7 × 10⁵ kg N y⁻¹ from land and 3.9 kg N ha⁻¹ y⁻¹ from the atmosphere to the Bay. Approximately 15% of the nitrogen load to the surface water and groundwater is due to the application of fertilizer in the watershed (Castro et al., 2003). In the northern section of the Bay, Baker et al. (2014) estimated that greater than 60% of the nitrogen load is from the Toms and Metedeconk rivers.

Sediment cores were collected on 29–30 July 2009 at four tidal marsh locations along the western shore of Barnegat Bay (Figure 1). Cores were collected during mid to low tide from interior locations away from any obvious disturbances (e.g., creek banks and ditches). At each site, two cores were obtained, one for chemical analysis and the other for stratigraphic descriptions of the marsh deposits. Push-piston cores approximately 1–1.5 m in length were collected using a tripod and pulley system. The cores were taken to the laboratory, extruded vertically, and sectioned in 2-cm intervals. Samples were stored in precleaned jars at −10°C.

In the laboratory, dry bulk density and loss on ignition (LOI) measurements were made for each core section to aid interpretations of radionuclide and nutrient data and to convert sediment depth (cm) to cumulative mass (g cm⁻²). Dry bulk density was calculated from gravimetric porosity using representative densities for pore fluid and organic and mineral solids. LOI was used to quantify the relative proportion of organic (combustible) and mineral (residual ash) solids in the core samples and was determined by combusting 4 g of sample powder in a muffle furnace at 550°C for 4 hours.

Sediment chronologies for the marsh sites were developed using down-core profiles of ²¹⁰Pb (t_1/2 = 22.3 y) and ¹³⁷Cs (t_1/2 = 30.1 y) activity concentration measured by gamma spectroscopy of the 46.5 and 661.7 keV photopeaks. Each of the 2-cm core samples were milled to fine powder, packed into a 60-mL plastic jar, and counted for 24−48 hours on Canberra Model 2020 low-energy germanium detectors. Excess ²¹⁰Pb (²¹⁰Pb_xs) was determined by subtracting the activity of a parent nuclide ²¹⁴Bi (609.3 keV) from the total activity (²¹⁰Pb_xs = ²¹⁰Pb_T − ²¹⁴Bi). Detector efficiencies were determined from counts of National Institute of Standards and Technology (NIST) Standard Reference Material 4357 (Inn et al., 2001). Because the counting geometry of the core samples and the NIST standard were identical, a self-absorption correction for ²¹⁰Pb was not necessary. Confidence limits reported with the radionuclide activities in this paper are the propagated one-sigma background, calibration, and counting errors. The computed minimum detectable activity (MDA) for ²¹⁰Pb and ¹³⁷Cs counted on the gamma detectors used for this study was 0.18 ± 0.01 dpm g⁻¹.

Profiles of ²¹⁰Pb_xs and ¹³⁷Cs activity were used to establish chronologies for the marsh sites. Marsh accretion rates were computed from ²¹⁰Pb_xs activity profiles following the simple constant flux–constant sedimentation (CFCS) model, which assumes that the depositional flux of ²¹⁰Pb and the sedimentation rate are constant (Appleby and Oldfield, 1992; Robbins, 1978). At steady state, the relationship between ²¹⁰Pb_xs activity and accretion rate is given by:

where, S is the accretion rate (cm y⁻¹), λ is the decay constant for ²¹⁰Pb (0.0311 y⁻¹), z is depth, A₀ is the activity of ²¹⁰Pb_xs at the marsh surface (dpm g⁻¹), and A_z is the activity of ²¹⁰Pb_xs (dpm g⁻¹) at depth z (cm). In this study, least-square regression of ln(²¹⁰Pb_xs) plotted vs. depth was used to determine the slope of the trend line, the second term on the right-hand side of Equation (1). Corresponding mass accumulation rates (MAR, g cm⁻² y⁻¹) were similarly computed by plotting ln(²¹⁰Pb_xs) against cumulative mass (CM). Rates computed in this manner eliminate the effect of variable porosity on the age-depth relationship. In the CFCS model, the age (t, y) of a sediment layer is related to MAR and CM as follows:

Confidence limits for ²¹⁰Pb-derived rates and ages reported in this paper were computed from the standard error of the slope of the regression line.

Mass accumulation and accretion rates were determined independently by ¹³⁷Cs chronology based on the depth of the activity peak and the assumption that the peak depth is concordant with 1964–the year of maximal atmospheric fallout from nuclear weapons testing (Ritchie and McHenry, 1990). The depth of the ¹³⁷Cs peak from the core top divided by the time between 1964 and the year of core collection gives the mean accretion rate, or mass accumulation rate when ¹³⁷Cs activity is plotted vs. cumulative mass. Confidence limits for the ¹³⁷Cs-based rates and ages were computed from the propagated radionuclide activity error. An advantage of ¹³⁷Cs chronology over ²¹⁰Pb is that it provides an absolute and potentially more accurate post-1964 chronology; however, ²¹⁰Pb chronology can provide a longer history of sedimentation, perhaps up to ∼100 years.

Cores were collected on the marsh surface (i.e. midmarsh area) during mid to low tide. Push-piston cores of approximately 1–1.5 m in length were retrieved by a tripod and pulley system. The cores were taken to the laboratory and sectioned into specific intervals (e.g., 2 cm). Samples were stored in precleaned jars at −10°C. Subsamples were thawed, dried at 60°C, and ground to a powder then stored in vials.

Total organic carbon and total nitrogen were measured using a CE Flash Elemental Analyzer following the guidelines in EPA 440.0 (Zimmerman, Keefe, and Bashe, 1997), manufacturer instructions, and the Academy's standard operating procedure. Samples were pretreated with acid to remove inorganic carbon. Total sediment phosphorus was determined using a dry oxidation method modified from Aspila et al. (1976) and Ruttenberg (1992). Solubilized inorganic phosphorus was measured with standard phosphate procedures using an Alpkem Rapid Flow Analyzer. Standard reference material (spinach leaves) and procedural blanks were analyzed periodically during this study. All concentrations are reported on a dry weight basis.

The stable isotopic composition of sediments was analyzed using a Finnigan Delta Plus coupled to an NA2500 Elemental Analyzer isotope-ratio mass spectrometer (EA-IRMS). Samples were run in duplicate or triplicate, with the results reported in the standard δ (‰) notation: δ_x = [(R_sample/R_standard) – 1] × 1000, where x is either ¹³C or ¹⁵N and R is either ¹³C/¹²C or ¹⁵N/¹⁴N. The δ¹⁵N standard was air (δ¹⁵N = 0), and for δ¹³C, the standard is the Vienna PeeDee Belemnite limestone that has been assigned a value of 0.0‰. Analytical accuracy was based on the standardization of the ultra–high purity N₂ and CO₂ used for continuous-flow IRMS with International Atomic Energy Agency (IAEA) N-1 and N-2 for nitrogen and IAEA sucrose for carbon, respectively. An in-house calibrated sediment standard was analyzed every 10th sample. Generally, precision based on replicate sample analysis was better than 0.2‰ for carbon and 0.6‰ for nitrogen.

The organic matter content of the core samples ranged from 10% to 70% by weight (Figure 2). This broad range reflects spatial and temporal variations in marsh plant productivity, organic matter decomposition, and burial of refractory organic matter with allochthonous mineral sediment among the coring sites. As is typical for Spartina marsh, organic matter content is highest within the living root zone of the sediment column, generally ∼10−20 cm below the surface. As in the general case, organic content varied inversely with dry bulk density as a reflection of differences between the densities of organic matter (1.2−1.5 g cm⁻³) and mineral particles (2.5−2.7 g cm⁻³).

Rates of sediment accumulation and accretion based on the ²¹⁰Pb_xs profiles ranged from 39 to 73 mg cm⁻² y⁻¹ and 0.16 to 0.29 cm y⁻¹, respectively (Table 1). Following Equation (2), the record of sedimentation provided by these profiles extended back to about ca. 1900. For all sites, ¹³⁷Cs activity increased up-core from the depth of first occurrence (15−25 cm) to a peak centered at 8−15 cm, above which activities decreased to levels close to MDA (Figure 2). Mass accumulation and accretion rates computed from the depth of the ¹³⁷Cs peak ranged from 34 to 84 mg cm⁻² y⁻¹ and 0.16 to 0.30 cm y⁻¹, respectively, consistent with rates derived by ²¹⁰Pb chronology.

Sediment organic carbon (SOC) concentrations for the four cores ranged from 3.2% to 33.7% on a dry weight basis, with an average of 15.7 ± 7.8% SOC (±1σ; Figure 3). The C to N ratio (atomic) of the marsh sediment ranged from 14 in the bottom sections of BB-4 to 41 in the upper section of BB-2 (Figure 3). In core BB-1, C:N at the surface was 17, increasing to approximately 22–28 (24 ± 4); similarly, in BB-4, C:N was approximately 16.8 in the surface 6 cm, increasing slightly to 22 then decreasing to between 14 and 20 (16 ± 2) below 18 cm. In the BB-2 site core, the sediment C:N value at the surface was 22, increasing to 42 at 10–12 cm then decreasing to between 22 and 32 (25 ± 3.6) below 12–14 cm. The core from BB-3 had a surface C:N value of 18, increasing with depth through the core to a maximum of 38 at 42–44 cm (Figure 3).

Total sediment nitrogen (TN) ranged from 0.23% to 1.4% N, with an average of 0.79 ± 0.34% (Figure 4), whereas total sediment phosphorus (TP) ranged from 0.019% to 0.39% P with an average of 0.071 ± 0.05% (Figure 5). TN concentrations were generally the highest in the upper 10 cm of the cores from sites BB-1, BB-2, and BB-4 (Figure 4). In these cores, sediment nitrogen contents below 20−30 cm were fairly constant, whereas at BB-3, sediment nitrogen content was variable throughout the depth. At site BB-3, in the discharge canal of the power plant, sediment nitrogen content ranged from 0.27% to 1.2%. In all cores, SOC distribution was also similar to TN. SOC concentrations were higher in the surface of cores from sites BB-1, BB-2, and BB-4 and decreased with depth. At sites BB-1 and BB-2, there were some slight variations at depth, whereas at BB-3, concentrations generally increased with depth from 15% at the surface to about 30% at 40 cm. Similar to C and N, TP decreased slightly with depth in all cores (Figure 5). The core from site BB-4 had a subsurface maximum of TP at 2–4 cm (0.39%), decreasing to near constant concentrations by 6–8 cm depth.

The carbon isotopic composition of the sediment ranged from −27‰ to −13‰ (δ¹³C average of −18.2‰ ± 4.5‰). Surface values were similar (ca. −18.5‰ to −16‰) and then increased slightly to approximately −14‰ to −16‰ within the root zone of cores from sites BB-1, BB-2, and BB-4. In the BB-3 core, the δ¹³C below the root zone decreased to −26.5 ± 0.6‰. Below the root zone (∼42–44 cm) in the BB-1 core, the δ¹³C became more negative (i.e. to −22‰), then increased to −15‰; beyond 44 cm depth, δ¹³C again decreased to −24‰. At BB-2 and BB-4, the δ¹³C remained fairly uniform with depth, with only slight variations (not shown here).

In cores from all locations, δ¹⁵N values ranged from −0.4‰ to 3.7‰. The δ¹⁵N in the sediment was higher in the surface cores, except in the BB-4 core, which had a higher signature belowground (Figure 6). In the BB-4 core, δ¹⁵N increased with depth from approximately 1.6‰ to 1.8‰ at the surface to 3.7‰ down-core. Except in the BB-4 core, δ¹⁵N increased from ∼30 cm belowground toward the surface. Nitrogen accumulation rates for the four sites ranged from 0.33 to 0.82 mg N cm⁻² y⁻¹ (Figure 7). Phosphorus accumulation rates ranged from 0.013 to 0.26 mg P cm⁻² y⁻¹. In all cores, rates were higher in the surface section of each core (0.04–0.16 mg P cm⁻² y⁻¹; Figure 8).

Excess ²¹⁰Pb activity in cores for all four sites decreased monotonically with depth, which indicates steady-state sediment accumulation and radioactive decay (Figure 2). The structure of the depth profiles shown in Figure 2 (decreasing bulk density with increasing LOI) indicate that the soil bulk density was influenced mostly by variation in the relative amounts of organic vs. mineral solids accumulation over time. Had postdepositional compaction been significant, soil bulk densities would have increased monotonically with depth, independent of organic content. Mass accumulation and accretion rates calculated from the ¹³⁷Cs activity peak indicated further that sediment accumulation rates at the marsh sites had not changed measurably after ca. 1954. Based on this result, and given that ²¹⁰Pb provides a longer history of sediment accumulation than ¹³⁷Cs, the ²¹⁰Pb CIC accretion rates were used to convert sediment depth to age to identify histories of nutrient burial. The rates of mass accumulation and accretion determined in this study are comparable to similarly determined rates reported for other locations within the Barnegat Bay system (Boyd, Sommerfield, and Elsey-Quirk, 2017; Unger et al., 2016), as well as tidal marshes elsewhere in the U.S. Mid-Atlantic region (Reed et al., 2008).

SOC, TN, and TP are higher than those found in tidal wetlands (salt and freshwater) in the Delaware River and Bay (Velinsky and Sommerfield, unpublished data; Velinsky, Charles, and Ashley, 2007; Velinsky et al., 2011). The C:N values increased slightly with depth (from 15 to 35 in the upper 20 cm) within the root zone of each site then decreased slightly (from 35 to 25) with depth (except for BB-3). The values suggest a mixture of Spartina-derived organic matter with some additional algal organic matter.

The carbon isotopic signatures identified from the cores indicate a Spartina-dominated system. Spartina is a C4 plant in which organic carbon produced would have carbon signatures generally around −12‰ to −14‰ (Long, Incoll, and Woolhouse, 1975). The δ¹³C signatures of organic carbon decrease with depth and were especially noted in BB-1 and to a lesser extent in BB-2, which may be due to decomposition of isotopically heavier polysaccharides that are preferentially degraded first (Benner, Fogel, and Sprague, 1991). Down-core changes in δ¹³C of organic matter in BB-3 may be due to habitat changes from nuclear facility construction and the building of the canal. Historical images from the 1930s (NJDEP, 2016) indicated that the current coring location (BB-3) was more of an upland site and not a tidal salt marsh. During the construction of the Oyster Creek Nuclear facility (mid-1960s), the creek was realigned, and the site became a wetland area now dominated by Spartina. Similar changes over time have been observed in the upper reaches of the tidal Murderkill River (Delaware Bay; Velinsky, Sommerfield, and Charles, 2010b) and northern San Francisco Bay tidal marshes (Byrne et al., 2001).

Changes in the δ¹⁵N of organic matter reflect potential changes in the source of N to the Bay. We could not identify the cause for the decrease in δ¹⁵N with depth, but it may be due to changes in the cycling and source of organic matter to these parts of Barnegat Bay. Cores obtained from farther downstream may help in this analysis. The stable isotopes of sediment nitrogen may give some indication as to changes in the source of nitrogen and biogeochemical cycling within the Bay. The higher δ¹⁵N is associated with higher inputs of N from urban sources and waste water (Ulseth and Hershey, 2005; Velinsky, unpublished data). Bratton, Coleman, and Seal (2003) identified similar δ¹⁵N trends in cores from the Chesapeake Bay. Similarly, Elliot and Brush (2006) showed a positive correlation between δ¹⁵N of tidal wetland sediments and nitrogen wastewater loadings over time. In the marshes of Woodbury Creek and Oldmans Creek, New Jersey, an increase in sediment δ¹⁵N was attributed to nitrogen discharge by wastewater treatment plants (Church et al., 2006). In the present study, an increase in δ¹⁵N (from <1‰ to between 3‰ and 10‰) at the surface started in the 1950s and was not directly related to the waste water from a facility that was built in the mid-1970s. Down-core variations of δ¹⁵N are somewhat difficult to interpret given the number of biogeochemical processes influencing nitrogen and variations of sources over time.

N and P accumulation rate (mg N or P cm⁻² y⁻¹) for each year were obtained by using the concentration of N and P in each interval. The missing intervals that were not analyzed for N and P were linearly interpolated. Even though the point source of nitrogen was redirected offshore, there was no sign of reduction in nitrogen accumulate rate after the redirection period (after 1975), perhaps because of continuous urban development in the Barnegat Bay watershed, indicating the contribution of other important sources (e.g., groundwater discharge and surface water runoff) to the Bay for the increase in nitrogen. In the upper bay core (i.e. BB-1), the rate doubled from the 1940s to the present (ca. 0.4–0.8 mg N cm⁻² y⁻¹). In the other three cores, the accumulation rate increased with time, but the increase was less pronounced.

At BB-1, the rate was variable over time (bimodal) with a subsurface peak around 1975, before the P ban in detergents and the redirection of point sources. The rate then decreased until the late 1990s but, afterward, exhibited a higher P accumulation rate. The BB-1 site, proximal to the Metedeconk and Toms rivers, would be exposed to higher nutrient loads that caused a high P accumulation rate after the 1990s. At BB-2, the rate of increase from the 1990s to present was similar to BB-1, whereas a higher accumulation rate from the 1990s to present was identified in BB-3 and BB-4.

Sediment profiles of N, P, and C reflect biogeochemical and physical processes for the diffusion for nutrients and sediment accretions. Most of the surface layer N and P in the Barnegat Bay varied because of the variation in river input. Many of these processes were substantially active in the upper sections of a marsh (e.g., root zones, 0–15 cm), and concentrations in this section may not reflect burial rate (generally concentrations of N and P are higher in the surface compared with depths >10–20 cm). As such, it is necessary to determine an average concentration of N and P in each core to account for diagenetic changes, as well as loading changes over time. For this, concentrations of N and P were multiplied by the dry sediment density (g cm⁻³) at each interval and then divided by the total mass of sediment that represents the past 60 years (i.e. 1950 to present). The average concentrations were then used along with the mass accumulation rate derived by ²¹⁰Pb chronology to provide an average accumulation rate for the past 60 years. The depth-integrated rates ranged from 4.5 to 6.0 g N m⁻² y⁻¹ (average = 5.2 ± 0.7 g N m⁻² y⁻¹) for nitrogen and from 0.33 to 1.0 g P m⁻² y⁻¹ (0.51 ± 0.33 g P m⁻² y⁻¹) for phosphorus and were slightly lower than those calculated for the surface section (3.3–7.6 g N m⁻² y⁻¹ and 0.43–1.6 g P m⁻² y⁻¹; Figures 7 and 8).

An estimation of the area of tidal coastal wetlands fringing Barnegat Bay is 26,000 acres (1.1 × 10⁸ m²; Lathrop and Haag, 2007; Rutgers School of Environmental and Biological Sciences, 2016). Using this area and the depth-integrated and core-top rates for N and P accumulation yield current burial rates (gross rates) of 5.48 and 6.45 × 10⁵ kg N y⁻¹, and 0.54 and 0.88 × 10⁵ kg P y⁻¹, respectively (Table 2). Using recent input estimates for nitrogen (Wieben and Baker, 2009) and phosphorus (BBNEP, 2005), coastal tidal marshes in the Bay can potentially sequester 79%–94% of the nitrogen and 54%–88% of the phosphorus (Table 2). Using the depth-integrated N and P rates over the past 50 years allows for the burial of 79 ± 11% and 54 ± 34% of the current input rate into the Bay. Similar burial rates and relative importance of marsh burial were shown by Greene (2005) for the Patuxent River (Maryland), by Velinsky, Sommerfield, and Charles, (2010b) for the Murderkill River (Delaware), and by Craft (2007) for marshes in coastal Georgia. These calculations show that marsh accumulation can sequester a majority of the P and N loads from the various sources (i.e. point and nonpoint sources). However, sediment recycling of N and P (Berner, 1980; Burdige, 2006) are not accounted for in these estimates and will modify and most likely reduce these estimates (i.e. Burial – Recycling = Net Burial). These estimates identify that the marshes and subtidal areas (Seitzinger, 1992) have a potential to trap both N and P before being exported to the Bay from the nontidal watershed and highlight the importance of ecosystem services that marshes provide (i.e. water filtration) and the potential cost of water treatment if marsh areas are reduced by either land development or sea-level rise.

Tidal wetland ecosystems are particularly vulnerable to rising mean sea level because they exist in a narrow vertical zone between mean tide levels and mean high water. Over the next 100 years, the rate of global mean sea-level (GMSL) rise is expected to accelerate from its current rate of 1.5−1.9 mm y⁻¹, which is based on the tide gauge network of GMSL (Church et al., 2013). As sea level rises or the marshland subsides, the seaward boundary of the marsh erodes, and the landward boundary extends inland when hydrogeomorphic conditions are suitable. However, in many coastal areas, topography and developed lands hinder marsh migration and thereby limit ecosystem services provided by tidal wetlands.

NOAA tide gauge records for Cape May, Atlantic City, and Sandy Hook, New Jersey, indicate local relative sea-level rise rates of 4.54 ± 0.55, 4.07 ± 0.16, and 4.05 ± 0.22 mm y⁻¹, respectively (NOAA, 2016). These rates are significantly higher than the GMSL trend because they reflect the combined effects of global and regional sea-level rise, as well as local land subsidence. Subsidence in Barnegat Bay area is reportedly due to collapse of the glacial forebulge and perhaps groundwater extraction (Davis, 1987; Engelhart et al., 2009). Using GPS data from 1995 to 2005, Snay et al. (2007) reported a local subsidence rate of 2.2 ± 0.7 mm y⁻¹ for the Sandy Hook tide gauge station. This suggests a local rate of absolute sea-level rise (4.05 − 2.2 = 1.85 mm y⁻¹) that is somewhat higher than the GMSL trend of 1.5−1.9 mm y⁻¹, perhaps because of regional ocean dynamics (Ezar et al., 2013). Subsidence rates for the other tide gauge stations are unknown, but they most likely fall in the range of 1−2 mm y⁻¹ on the basis of data presented by Snay et al. (2007). Indeed, for the New Jersey coast, Engelhart et al. (2009) reported subsidence-corrected rates of sea-level rise of 1.3 ± 0.2 and 1.4 ± 0.7 mm y⁻¹, averaged over the past 4000 years. These rates are somewhat lower than those suggested by GPS data for the coast of New Jersey (Snay et al., 2007), but this may be related to differences in averaging periods and temporal variation in subsidence rates.

Marsh accretion rates determined by ¹³⁷Cs and ²¹⁰Pb chronology in this study ranged from 1.6 ± 0.2 to 3.0 ± 0.4 mm y⁻¹. Assuming the marsh landscape as a whole is subsiding at 2.2 mm y⁻¹ (Snay et al., 2007), then the corresponding change in marsh surface elevation (accretion minus subsidence) in the study area ranges from −0.6 mm y⁻¹ (BB-2) to 0.8 mm y⁻¹ (BB-4). These rates fall below the aforementioned tide gauge rates of sea-level rise (1.85−2.34 mm y⁻¹) upon removing the effect of land subsidence (2.2 mm y⁻¹). Hence, the apparent accretionary deficit renders these marshes vulnerable to inundation, particularly should the rate of sea-level rise accelerate in the future. Here it should be noted that subsidence is highly variable on small spatial scales, so the 2.2 mm y⁻¹ value used here is speculative and perhaps unrepresentative of actual conditions. Nonetheless, even if subsidence is as low as 1 mm y⁻¹, these marshes are accreting at below local rates of relative sea-level rise. An ongoing effort to identify local subsidence rates in the marshlands should provide data needed to place the accretion rates in the context of sea-level rise.

This study involved the chemical analyses of sediments subsampled from four sediment cores taken from tidal marshes in Barnegat Bay, New Jersey. Four cores were collected along the western coast of the Bay and dated using ²¹⁰Pb and ¹³⁷Cs chronology. All cores provided temporal coverage (>50–80 y) for detailed chemical analyses. Sediment N accumulation rates have increased slightly from the early 1960s to the present. The trend in P accumulation was variable between sites. Only in the upper bay core was there a subsurface maximum in the mid-1970s at the time of the Clean Water Act ban in P detergents (Church et al., 2006). All cores showed higher accumulation rates nearer the surface. Mass balance calculations show that the remaining wetlands in the Bay can sequester, via burial, approximately 79% N and 54% P of the incoming watershed load. This highlights the ecosystem services that tidal wetlands provide to the Bay and coastal waters. Marsh accretion rates ranged from 0.16 to 0.30 cm y⁻¹. Comparison of marsh accretion rates and tide gauge records of local relative sea-level rise suggests that Barnegat Bay marshes are vulnerable to inundation should the rate of sea-level rise accelerate in future.

This work was funded by New Jersey Department of Environmental Protection and Patrick Center for Environmental Research, Academy of Natural Sciences at Drexel University. The authors thank Roger Thomas, Jared Halonen, Paul Kiry, and Paula Zelanko for field and laboratory assistance and data interpretation and Izabela Wotjenko (U.S. Environmental Protection Agency) for her valuable support. The authors are grateful to two anonymous reviewers for their suggestions to improve the manuscript.