Dataset: Alkalinity, salinity, bivalve biomass, streamflow, and submerged aquatic vegetation.

The data are subsets from other data sets: the Chesapeake Bay Program's (CBP's) Water Quality Database, United States Geological Survey (USGS) monitoring data and Virginia Institute of Marine Science (VIMS) long-term aerial imagery monitoring dataset for the Chesapeake Bay (see related datasets). Here we provide a description of the subsetting and any quality control.

Detailed description of estuarine water chemistry data and processing

All alkalinity measurements made in tidal waters of the Chesapeake Bay were downloaded from the Chesapeake Bay Program’s Water Quality Database (https://www.chesapeakebay.net/what/downloads/cbp_water_quality_database_1984_present) in March, 2019. A total of 26,504 alkalinity measurements across 95 stations were identified. Total alkalinity was measured by titrating a water sample to a pH of 4.5 (Chesapeake Bay Program, 1996), a method that has a precision of 0.02 mol m–3 (Strickland & Parsons, 1972). Reported as mg CaCO3 L–1, alkalinity was converted to mol m–3 by assuming 100 g of CaCO3 corresponds to 2 moles of alkalinity. (Note that, when discussing alkalinity quantitatively, dimensions of concentration are implied unless otherwise noted.) Seven data points were removed from the data set, six being identified as extreme outliers in this study and one flagged by the Chesapeake Bay Program. Outliers were defined based on the entire dataset as values more than six interquartile ranges below or above the 25th and 75th percentiles, respectively. 

Salinity information is helpful in the analysis of alkalinity data and thus we wanted to pair every alkalinity measurement with either a salinity measurement or, in some cases, a salinity value of zero if salinity was not measured but we were confident that the water was fresh. Before pairing alkalinity and salinity measurements, any replicate measurements reported for a given time and depth were averaged (1,140 measurements for alkalinity and 129 for salinity). A large fraction, 26%, of the alkalinity measurements did not have corresponding salinity measurements. Nine of the stations without salinity measurements had only one alkalinity measurement and were thus deleted. For the remaining stations (which included stations in tidal fresh water and saline water), an alkalinity measurement was removed if there was no simultaneous salinity measurement and if the mean salinity at the station was greater than 0.5. We removed the data because we did not want to include in the analysis an alkalinity measurement whose corresponding salinity was not known. For the remaining alkalinity data without corresponding salinity data, the salinity was set to zero. Finally, we eliminated a total of 6 stations along the eastern and western shores of the mainstem bay, which contained a total of only 12 observations. This processing reduced the alkalinity data set to 25,289 measurements across 80 stations (“whole processed dataset,” hereafter). 

For certain analyses, a subset of the whole processed data set (“reduced data set,” hereafter) was created so that robust mean annual cycles and multi-year averages could be computed over a common time period in each tidal tributary. The Chester River Estuary was not included because it contains only one station and does not have corresponding nontidal alkalinity data. Stations with good seasonal coverage (at least monthly resolution) over at least 3 years were selected and data for those stations outside of a common period of at least two years were removed. Data outside of the main channels of the tidal tributaries were also excluded, notably stations in tributaries of the Potomac River Estuary, such as the Anacostia River, Piscataway Creek, and Mattawoman Creek. The last five years of the Potomac River Estuary record were also removed because no riverine alkalinity data were available for that period. Finally, to facilitate comparisons among the tidal tributaries, the data set was limited to the upper 5 m. Across the seven tidal tributaries, the number of stations was reduced from 80 to 53 and the number of alkalinity observations from 25,289 to 11,007. An additional quality control step was included to prevent potential biases in the mean annual cycles. Alkalinity outliers in this step were defined as above for the full data set but using more stringent criteria: on a station-by-station basis, rather than all stations together, and using a threshold of 2, rather than 6, interquartile ranges. The result was the removal of another 22 outliers. Mean annual cycles were created from this data set by first averaging all data (in the upper 5 m) within a given month of a given year. From 1 to 13% of the station-months had gaps, which were filled using linear interpolation. The resulting depth-averaged monthly gridded dataset contained 5,788 data points. Then all Januarys were averaged, all Februarys were averaged, etc. The multi-year average was computed by averaging the months of the mean annual cycle.

To better understand sources and sinks of alkalinity, we quantified the impact of nitrogen cycling by exploiting measurements of nitrate, nitrite, and ammonium from the same source that the alkalinity measurements were taken from; measurement methods are described in Chesapeake Bay Program (2012). Any of these nitrogen data flagged by the Chesapeake Bay Program were removed, except for data below the detection limit, which were set to zero. 81.4 and 99.7% of the alkalinity data in the whole data set and the reduced data set, respectively, had paired estimates of A_N.   

Detailed description of riverine water chemistry data and processing:

Riverine alkalinity and streamflow data were derived from seven gauging stations of the United States Geological Survey (USGS, https://nwis.waterdata.usgs.gov/nwis, accessed in March 2019). The gauges for the Potomac and Rappahannock Rivers are located just a few km above the fall line, those for the Susquehanna and James Rivers are well above the fall line, and those for the Patuxent, Mattaponi, and Pamunkey Rivers are well below the fall line. Insufficient data were available for the Chester River to conduct an analysis of its flow and alkalinity. The Weighted Regressions on Time, Discharge, and Season (WRTDS) model (Hirsch et al. 2010) was used to estimate daily concentrations and fluxes of alkalinity. The daily fluxes were summed for each month to get monthly fluxes and then an effective monthly mean alkalinity concentration was determined by dividing the monthly alkalinity flux by the monthly streamflow. The effective monthly mean alkalinity concentration is what the estuary experiences on average over the course of a month because the impact of the river is felt more during high-streamflow events than during low-streamflow events. An effective long-term mean alkalinity concentration was also computed by dividing the long-term average flux by the long-term average streamflow. 

The riverine data cover the time periods of the estuarine data except for the years 2014–2018 for the Potomac River and the last three months of 1999 for the Rappahannock, Mattaponi, Pamunkey, and James Rivers. For the latter four rivers, the last three months of 1999 were filled in using a least squares model of monthly alkalinity as a power-law function of monthly streamflow based on data from January 1996 to September 1999. 

The gauges do not capture all nontidal inputs of freshwater and alkalinity to the tidal tributaries. We used online USGS estimates of areas of hydrologic units and the USGS StreamStats application to estimate the area of the watershed draining to a particular station in a tidal tributary and what fraction of that area is monitored by the gauge available in that watershed. The gauges capture between 40 and 87% of the watershed draining to the most downstream station in each tidal tributary. The total area gauged is 76% of the Chesapeake Bay watershed.

To quantify the influence of nitrogen cycling on riverine alkalinity, estimates of the monthly mean nitrate + nitrite concentration at the seven gauging stations were acquired from USGS’s River Input Monitoring (RIM) stations for the Chesapeake Bay (Moyer & Blomquist, 2018). Estimates of nitrate + nitrite concentration are based on the WRTDS model of Hirsch et al. (2010) and cover the period 1984 to 2017. Unfortunately, monthly mean estimates of ammonium concentration were not available. To assess the contribution of ammonium, we analyzed limited ammonium measurements available from the Chesapeake Bay Program  (http://data.chesapeakebay.net/WaterQuality, accessed March 2019) at the gauging stations for the period 2012 to 2017, when sufficient ammonium data were available for comparison with the RIM nitrate + nitrite product. The nontidal data were downloaded from the same source that the tidal alkalinity and nitrogen measurements were taken from. Quality control was applied to the ammonium data by removing flagged data, except for those that had detection limit flags. Data flagged as being below the detection limit were set to zero. Some data were also flagged as being above the detection limit, but our viewing of these data made it clear to us that these data were, in fact, below the detection limit (i.e., a transcription error was made). We therefore also set to zero data flagged as being above the detection limit.

Processing and analysis of bivalve biomass and SAV coverage

The bivalve biomass data are Chesapeake Bay Program measurements of ash-free dry mass, which is an estimate of organic (i.e., nonmineral) biomass per unit area. Details of the methodology are given in Llansó and Zaveta (2017, and references therein) and are summarized here. The Chesapeake Bay Program uses two sampling strategies: fixed stations that are revisited typically once each year and random stations that differ in location from year to year. Since 2009, both fixed and random station sampling has occurred during the summer months, although additional months were sampled historically. Two fixed stations (36 and 40) are located within the boundary of the Potomac box model. Station 36 is located in tidal fresh waters and station 40 in oligohaline waters. The number of months that sampling occurred changed over time but there was at least one sampling effort made between August and October of every year the box model was run (1986–2013), except for 1992 for station 36 and 1989 for station 40. For each month that a station was visited, triplicate benthic grabs were collected; these were averaged for a given month. If multiple months in the August–October period were sampled in a given year, which occurred seven times at each station, the monthly means were averaged to produce a single value for the August–October period.

Estimates of bivalve calcification rate were made from measurements of ash-free dry mass per unit area (Bdry) following the procedure of Chauvaud et al. (2003), which has three steps: (1) conversion of Bdry to organic carbon biomass per unit area (B); (2) conversion of B to minimum potential secondary production per unit area (Pmin), an estimate of the organic matter production rate; and (3) conversion of Pmin to the calcification rate per unit area (Cmin). For the first step, a conversion factor of 0.41 organic carbon mass per ash-free dry mass was used. For the second step, a specific growth rate of 4.45 yr–1 was used, which is typical for the dominant bivalve in low-salinity waters of the Potomac River Estuary, Corbicula fluminea. This specific growth rate corresponds to a turnover time of 82 days and is derived from McMahon (2002), who gives a range in Corbicula turnover times of 73–91 days. For the third step, a ratio of 15 g CaCO3 per g C (ratio of shell production to tissue production) was used. This value is probably the most uncertain and is based on unpublished data cited by Chauvaud et al. (2003). In summary, there are three conversion factors needed to compute calcification rate from ash-free dry mass. The resulting equation for calcification rate is: Cmin = 27.4 g CaCO3 (g dry mass)–1 yr–1 Bdry.

To determine the representativeness of the two fixed benthic monitoring stations, a comparison of Bdry was made with the stations that were sampled randomly. The random sampling strategy began in 1995 and the comparison was limited to 1995–2013. Random sampling during this period ranged from August to October. Polygons corresponding to the tidal fresh (Box 1) and oligohaline (Box 2) regions were constructed using ArcGIS (v10.4, ESRI) and all random stations that occurred in a given polygon were identified and extracted. Annual averages were calculated for each box and year from the available stations, with the exception of 1998 and 2003 for Box 1 and 1999 and 2006 for Box 2, years for which no stations fell within the respective polygons. For years in which stations did fall within each box, the number of stations ranged from 1 to 8 per year, averaging 2 per year for Box 1 and 3 per year for Box 2.

Areal coverage of SAV was estimated for Boxes 1 and 2. SAV data were from the Virginia Institute of Marine Science long-term aerial imagery monitoring dataset for the Chesapeake Bay (Orth et al., 2018). Measurements are made at a given location annually, typically between June and October, to capture peak conditions (which vary with species). The total area of SAV, regardless of SAV density classification, was used. The area of SAV in each region was obtained for each year from 1984 to 2017 except for 1988 and is projected at a scale of 1:24,000. Using ESRI ArcGIS (ESRI, 2011), we first clipped SAV coverages to the management segmentation scheme boundaries used by the Chesapeake Bay Program that correspond to the tidal fresh and oligohaline regions of the Potomac River Estuary (POTTF and POTOH) and then computed corrected polygon areas to calculate total area of SAV in each box.