Age variability and decadal time-averaging in oyster reef death assemblages

Decades of work on molluscan death assemblages (DAs) have successfully documented temporal changes in community composition or species attributes from direct assessments of the remains themselves (e.g., Kowalewski et al., 2000; Kidwell, 2007; Dietl and Durham, 2016; Albano et al., 2021) or from proxy information derived from them (e.g., Gillikin et al., 2019). Despite the promise of such geohistorical records for conservation paleobiology, examples of their use by resource managers are still uncommon (Groff et al., 2023). One reason is the difficulty of putting DA data in temporal context. Geochronological analyses (e.g., radiocarbon dating) are expensive and difficult to interpret, leading many conservation paleobiological studies to work around age-related uncertainties by citing general assumptions and/or studies from similar depositional settings (e.g., Dietl and Durham, 2016).

However, assemblage- or specimen-level chronological control is often required to meaningfully compare DA data with the annual or subannual real-time monitoring data typically used for resource management. This was the case for the Historical Oyster Body Size (HOBS) project in Florida, southeastern United States—codeveloped by the Florida Department of Environmental Protection (FDEP) Office of Resilience and Coastal Protection (ORCP) and the Paleontological Research Institution (PRI; Dietl et al., 2023)—which aimed to use oyster reef DA samples to supplement real-time monitoring data on oyster body sizes for ORCP’s Statewide Ecosystem Assessment of Coastal and Aquatic Resources (SEACAR) project (www.floridadep.gov/SEACAR).

The aquatic preserves managed by ORCP were established between 1966 and 2020 to be maintained “in an essentially natural or existing condition” (Florida Administrative Code R.18-20.001[2]; Florida Department of State, 1997). Thus, management of each preserve is often focused on its relative condition since establishment, meaning the ultimate utility of the DA approach for SEACAR would be influenced by the specific age and time-averaging properties of the oyster reef DAs. We hypothesized that oyster reef structure might limit postburial stratigraphic mixing enough such that samples from the DAs could yield data at a high enough temporal resolution to be integrated with real-time monitoring data from living oyster populations. To test this assumption and develop an understanding of both oyster reef taphonomy and the potential utility of DA data for FDEP, we produced a geochronological data set to quantify the absolute ages and temporal resolutions of oyster reef DAs from around the state.

Here, we describe this investigation and show that oyster reef DAs preserve reliably recent and high-resolution stratigraphic records relative to most other molluscan DAs documented to date, suggesting these records are often appropriate for decadal-scale conservation paleobiological investigations. We also highlight the geographic variability in our data set and its implications for the importance of location-specific geochronological information for increasing the salience of paleoecological data for the resource management community.

In order to build a geochronological data set to evaluate the utility of oyster DA samples for documenting trends over recent decades, we randomly selected 630 Crassostrea virginica left-valve specimens from oyster DA samples representing two stratigraphic intervals (15–25 cm and 25–35 cm) collected from up to three sample holes positioned across the densest living portion of each of 31 natural, intertidal oyster reefs in 11 locations around Florida (see the Supplemental Material¹), i.e., between 2 and 7 specimens from each DA sample (Fig. 1). The selected specimens were dated by radiocarbon analysis of powdered carbonate targets (Bush et al., 2013; Hua et al., 2019)—a less expensive method with lower precision than the standard analysis of graphite targets, but one that yields similar ages (Bright et al., 2021)—to achieve a higher sample size (see the Supplemental Material for details on specimen selection for radiocarbon analysis as well as a sample size validation using 80 additional randomly selected specimens from four of the DA samples). Specimens were prepared at Northern Arizona University (NAU; Flagstaff, Arizona, USA) and analyzed at either the W.M. Keck Carbon Cycle Accelerator Mass Spectrometry facility at the University of California, Irvine, or NAU’s own Arizona Climate and Ecosystems (ACE) Isotope Laboratory. Local corrections for the hardwater effect (e.g., Spennemann and Head, 1998) and/or estuarine influences (e.g., Ulm et al., 2009), in terms of dead carbon contribution, were developed using additional radiocarbon analyses of two live-collected oyster specimens from each sampling area (see the Supplemental Material).

Age calibration was performed using OxCal v4.4 software (Bronk Ramsey, 2009) and the Marine20 calibration curve (Heaton et al., 2020) with a constant regional marine reservoir correction, ΔR = −134 ± 26 yr, which is equivalent to 5 ± 32 yr (Kowalewski et al., 2018) relative to Marine13 (Reimer et al., 2013), extended to 2022 using a regional marine bomb radiocarbon curve based on our data as well as 665 other radiocarbon results from the Gulf of Mexico, western Atlantic Ocean, and Caribbean Sea from 24 additional studies (see the Supplemental Material). Following Kowalewski et al. (2018), we used empirical posterior distributions of age probabilities for the specimens in each DA sample to generate estimates of (1) DA sample ages (we use the terms “specimen age” and “sample age” to refer to radiocarbon results for an individual oyster shell and all oyster shells from a given DA sample, respectively), and (2) time-averaging. Due to recently published concerns about the corrected posterior age estimate (CPE; sensu Kowalewski et al., 2018; also known as residual time-averaging in some studies), however, we used the interquartile range (IQR) of the average sample age probability distribution, with the quartiles weighted by the age probabilities—the total age variability (IQR_TAV)—alone to estimate time-averaging instead of the IQR_TAV and CPE (Ritter et al., 2023; see the Supplemental Material).

Finally, to compare the contributions of location and burial depth to overall variation in DA sample median age and IQR_TAV, we fit a hierarchical Bayesian model to the data for each burial depth as well as the burial depth difference for each DA sample hole (see the Supplemental Material). All data analyses were conducted using R statistical software v4.3.0 (R Core Team, 2023) and RStudio (RStudio Team, 2023).

The radiocarbon results indicated that oyster reef DAs are high-resolution archives with abundant shells from the recent past and minimal time-averaging in comparison to other molluscan DAs. Among the 126 dated oyster DA samples, median calibrated ages ranged from 1567 to 2012 CE, but 91% were post-1950 (Fig. 2), and 6.4% of the DA samples had subdecadal-scale IQR_TAV (0–10 yr), 72.8% had decadal-scale IQR_TAV (11–100 yr), and 20.8% had centennial-scale IQR_TAV (101–1000 yr) (Fig. 3; see Appendix S1 in the Supplemental Material for DA sample-level results). Moreover, collocated samples from different burial depths showed the expected temporal order (i.e., deeper = older) in most cases: Out of the 53 sample holes for which both depth intervals were processed and dated, 12 had median DA sample ages for the 15–25 cm burial depth that were older than those of material from the 25–35 cm burial depth, and five of those cases were from a single locality (Lone Cabbage; Fig. 2). The results also showed that the age and time-averaging of a given burial depth can vary substantially over small spatial scales (i.e., both intrareef and interreef assemblage variation; Fig. 2). In fact, the modeled standard deviations (SDs) for spatial variability in median age and IQR_TAV (e.g., DA sample-hole-level median SDs were 12.8 and 20.9 yr for median age and IQR_TAV, respectively, for the 15–25-cm-depth samples) were of similar magnitude to those for the difference between burial depths (e.g., DA sample-hole-level median depth difference SDs were 29.7 and 37.0 yr for median age and IQR_TAV, respectively; Table S5; Figs. S5–S10).

To our knowledge, this is the largest study of age-depth relationships and the first study to use radiocarbon to document time-averaging in oyster reef DAs. We found that, relative to other molluscan DAs, the oyster DA samples were younger, were less time-averaged, and had less spatial variability in both calibrated age and time-averaging estimates (Flessa et al., 1993; Meldahl et al., 1997; Kowalewski et al., 1998, 2018; Kosnik et al., 2009, 2015; Krause et al., 2010; Dexter et al., 2014; Dominguez et al., 2016; Ritter et al., 2017; Tomašových et al., 2019; Albano et al., 2020; see additional studies summarized by Kidwell [2013, their table 1]; but see also Tomašových et al. [2018] for an example of a non-reef DA with decadal-scale resolution). This result agrees well with that of a preliminary investigation of time-averaging on two southwest Florida oyster reefs by Lindland et al. (2001) that used amino acid racemization geochronology.

Among the few recent studies that (1) focused on mollusks, (2) estimated time-averaging and sample ages in similar ways to our study, and (3) reported their unsummarized, sample-level results, the C. virginica DA samples typically had younger median ages by ~100 yr, and over half of our samples also had lower IQR_TAV, some by an order of magnitude or more (Fig. 3). For instance, Dominguez et al. (2016) sampled the upper 20 cm of sediment (medium–fine sand, <2% mud) at six sites with ~9 m water depth in Sydney Harbor, Australia, and found time-averaging (IQR_TAV) of ~84–>2000 yr in DAs of the bivalve Fulvia tenuicostata, but the median ages of the samples were ~150 yr. In contrast, IQR_TAV across all of the C. virginica DA samples in our study ranged from 6 to 532 yr with a median of ~25 yr. The medians of the median calibrated ages across all of the DA samples from 15–25 cm and 25–35 cm burial depths were 23 yr and 29 yr, respectively. The SDs for both median age and IQR_TAV among the locations sampled by Dominguez et al. (2016) were higher than the respective modeled locality-level SDs for the oyster reefs we sampled, despite the much greater geographic area covered by our study (Figs. S5–S10). A similar pattern is evident for DAs from other depositional settings and locations (e.g., southern Brazilian shelf—Ritter et al., 2017; subtidal sand flat, Fernandez Bay, San Salvador Island, Bahamas—Kowalewski et al., 2018; Fig. 3).

One exception to this pattern is the study by Tomašových et al. (2018), which found comparable age and time-averaging estimates to ours in Corbula gibba DAs from cores of the Po and Isonzo prodeltas, northern Adriatic Sea (Fig. 3). However, the authors stated that the two deltas have some of the highest sedimentation rates in the northern Adriatic Sea, and median ages and time-averaging estimates for C. gibba DAs from the eastern Gulf of Trieste—across the gulf from the Isonzo River and characterized by low sedimentation rates—were older and more time-averaged than the prodelta samples by nearly two orders of magnitude (Tomašových et al., 2019). In contrast to these large differences in age and time-averaging of C. gibba DAs between depositional settings, decadal-scale resolution appears to be a common feature of DAs from intertidal C. virginica reefs in multiple estuaries across Florida.

Overall, our results suggest that oyster reefs have a relatively high shell burial rate and less stratigraphic mixing relative to nonreef molluscan DAs, consistent with the hypothesis that the physical structure of oyster reefs reduces the susceptibility of DAs to some taphonomic processes. Despite their higher temporal resolution than other types of molluscan DAs in most cases, there is still considerable variation in the oyster DA median ages and IQR_TAV values (Fig. 2), precluding useful regional or statewide generalizations of age versus burial depth relationships or scales of time-averaging (see the Supplemental Material for an example).

This variability highlights a need for additional work to further refine the spatial and temporal specificity of dead C reference information for calibrating the fraction modern carbon (F¹⁴C) of estuarine carbonates to improve their accuracy and precision as much as possible (see Supplemental Material for further discussion of this point). It also illustrates why specific geochronological information will be important for many conservation paleobiological contributions to oyster management. Exactly how necessary they are for any given project will depend on the questions investigated, but trends in many indicators of oyster population condition, such as live-oyster size-frequency, are typically tracked at annual or subannual intervals by oyster monitoring programs. To integrate measurements from DA samples with such high-resolution records for trend analyses, it would likely be necessary to know, for instance, whether the median calibrated age and IQR_TAV of a DA sample are 2011 and 53 yr, respectively, or 1979 and 16.5 yr—as was the case for two of the DA samples at 15–25 cm burial depth from our Big Hickory locality.

Once these data are obtained, comparisons between the DA data and monitoring data that were impractical without them can become feasible—such as the HOBS project’s focus on integrating DA and real-time oyster size data into a single trend analysis, accounting for uncertainty in both oyster size and sample ages—instead of only focusing on more general “before/after” comparisons (e.g., Dietl and Durham, 2016). Further, our study demonstrated that most of the DA samples from oyster reefs represent a relevant time period for ORCP management (i.e., late 1960s to mid-2000s) and can yield decadal-scale (and sometimes subdecadal-scale) retrospective information from ORCP-managed areas where no long-term contemporaneous oyster monitoring took place.

Lastly, the apparently limited stratigraphic mobility of shells in recent oyster DAs suggests the intriguing possibility that the degree of time-averaging in an in situ fossil oyster reef (bioherm) is not dramatically greater than that in the DA of a living oyster reef. In this case, bioherms might preserve decadal-scale records from time periods when information at such a fine temporal resolution is exceptionally rare, making them suitable records for otherwise impossible studies of short-term ecological processes in the deep past (e.g., Kowalewski et al., 1998; Kidwell and Tomašových, 2013).

We thank Jordon Bright (Amino Acid Geochronology Laboratory, Northern Arizona University) and the staffs of the Keck Accelerator Mass Spectrometer Facility (University of California, Irvine) and the ACE Isotope Laboratory (Northern Arizona University) for ¹⁴C analyses, and Matthew Kosnik (Macquarie University, Australia) for help in adapting the Kowalewski et al. (2018) R script. We also thank colleagues and volunteers from the Florida Department of Environmental Protection (FDEP), the Paleontological Research Institution (PRI), and the Florida Department of Agriculture and Consumer Services (FDACS) who helped with field or laboratory work (*FDEP, ^†PRI, ^‡FDACS, ^§Pelican Island National Wildlife Refuge). Staff: M. Anderson*, P. Benjasirichai^‡, E. Bourque*, M. Brown*, C. Brunk*, C. Clark^‡, S. Cofone^‡, R. Cray*, E. Dark*, M. DeHaven^‡, N. Dix*, S. Erickson*, J. Fleiger^‡, J. Garwood*, T. Green*, B. Hamill*, K. Harshaw^‡, D. Hersl*, T. Jones*, K. Lang*, P. Marcum*, M. McMurray*, B. Mowbray*, R. Noyes*, J. Pier^†, R. Prado*, M. Pruden^†. Volunteers: B. Alexander*, A. Bishop*, S. Gavirneni^†, C. Hormuth*, A. McNeil^†, M. Melekos^†, R. Mondazzi*, D. Philipp*, E. Prest^†, M. Raymond*, M. Schilling*, Z. Siper^†, B. Skoblick^†, A. Thorsness*, D. Thorsness*, J. Valentine^§. We also gratefully acknowledge K. Flessa, M. Savarese, and an anonymous reviewer, whose comments helped to improve the manuscript, and funding from the National Oceanic and Atmospheric Administration (NOAA) Office for Coastal Management under the Coastal Zone Management Act of 1972, as amended, to the Florida Coastal Management Program (NOAA awards NA18NOS4190080 and NA19NOS4190064). The views, statements, findings, conclusions, and recommendations expressed herein are those of the author(s) and do not necessarily reflect the views of the State of Florida, NOAA, or any of their subagencies.