Using sedimentological priors to improve 14 C calibration of bioturbated sediment archives

Radiocarbon (C) dating is often carried out upon multi-specimen samples sourced from bioturbated sediment archives, such as deep-sea sediment. These samples are inherently heterogeneous in age, but current C calibration techniques applied to such age heterogenous samples were originally developed for age homogeneous material. A lack of information about age heterogeneity leads to a systematic underestimation of a sample's true age range, as well as the possible generation of significant age-depth artefacts during periods of highly dynamic ΔC. Here, a calibration protocol is described that allows for the application of sedimentological priors describing sediment accumulation rate, bioturbation depth and temporally dynamic species abundance. This Bayesian approach produces a credible calibrated age distribution associated with a particular laboratory C determination and its associated sedimentological priors, resulting in an improved calibration, especially in the case of low sediment accumulation rates typical of deep-sea sediment. A time-optimised computer script (biocal) for the new calibration protocol is also presented, thus allowing for rapid and automated application of the new calibration protocol using sedimentological priors. This new calibration protocol can be applied within existing age-depth modelling software packages to produce more accurate geochronologies for bioturbated sediment archives. 1.0 Introduction Radiocarbon (C) analysis is routinely used to determine the age of marine sediment archives, and has been fundamental in increasing understanding the spatio-temporal development of global ocean and geochemical parameters during the last glacial and the Holocene. However, due to C being a very rare radioisotope in the environment (approximately one in one trillion carbon dioxide molecules in the atmosphere is CO2), it is more difficult to measure than more common, stable carbon isotopes. From a practical standpoint, this rarity results in a requirement of relatively large 5

This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30.
Feedback is welcome, please contact the author.
microfossils within a sediment archive results in a systematically incorrect calibration, which can lead to the production of age-depth artefacts . These artefacts are further amplified for sediment intervals coinciding with highly dynamic Δ 14 C and/or Δ 14 C plateaus (ibid.).
Such age-depth artefacts could potentially be misidentified as true age-depth features in the sediment, and incorrectly attributed to, e.g., environmental and/or climatological processes.

Method
The calibration process presented here involves using quantitative priors for sediment accumulation rate, bioturbation depth and temporal changes in species abundance to estimate the credible age distribution for a sample that would result in the observed 14 C activity determination derived from a given sample. In addition to the prescribed priors, this process must also take into account all 14 C uncertainties, i.e uncertainties pertaining the laboratory determination, past Δ 14 C (the calibration curve) and reservoir effect. To avoid ambiguity, throughout this text the use of the term "age" refers exclusively to true/calibrated age, while 14 C activity is always referred to as 14 C activity, i.e. not as " 14 C age".

Establishing a prior distribution for calendar age
In order to calibrate 14 C activity measurements carried out upon heterogeneous samples retrieved from bioturbated sediment, the following sedimentological priors are defined: s = estimated sediment accumulation rate (SAR), in cm yr -1 m = bioturbation (mixing) depth (BD), in cm. k = the fraction of the analysed microfossils that are fragmented (a value between 0 and 1) a = time series of abundance of the analysed species relative to itself (values between 0 and 1) Both SAR and BD are considered here as a single value, i.e. not as a time series of temporally variable values. These parameters are kept static foremost to reduce computation time, and also because temporal changes in, e.g., SAR (the relationship between mean age and depth) are not known when an age-depth chronology has yet to be developed. In short, applying detailed information about temporal changes in SAR when the age of the sediment is not yet known would constitute circular thinking. This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30. Feedback is welcome, please contact the author.
The sedimentological priors can be used to construct a prior distribution of relative age for the sample being calibrated. Following (Berger and Heath 1968), the age distribution for a given depth of fully bioturbated sediment core can be represented by an exponential probability distribution, which can be considered the basis of the prior probability distribution for our sample's calibrated age: where r is the relative age (starting at 1 yr) within P prior . The low-probability long tail of an exponential probability function continues to infinity, which obviously cannot be stored in computer memory. The prior distribution is therefore limited to the age equivalent value of five bioturbation depths, i.e. a relative age of r limit = 5m / s, which is rounded to the nearest whole year.
When picking microfossils for 14 C analysis, palaeoceanographers generally prefer to pick whole specimens. The fragmented and/or dissolved microfossils that are not picked have been exposed to more bioturbation samples, and as such represent the oldest fraction of the sample (Rubin and Suess 1955;Ericson et al. 1956;Emiliani and Milliman 1966;Barker et al. 2007). Information regarding this effect should be incorporated into the prior distribution. The fraction of fragmented microfossils (k) can be related to the cumulative expression of Eq. 1 as follows: Eq. 2 can be solved to attain r(k), the threshold age for fragmented foraminifera: This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30. Feedback is welcome, please contact the author.
process. When r(k) ≥ r limit , r(k) is approximated to r limit . All discrete probability values in p prior are subsequently normalised such that they sum to 1.

Establishing a distribution for 14 C activity
The calibration process carried out must incorporate the full uncertainty regarding 14 C activity, which includes uncertainties regarding the laboratory 14 C activity determination, the calibration curve 14 C activity, and the 14 C activity depletion as a result of the reservoir effect. These are expressed here as follows: A det = The laboratory 14 C activity determination of the sample (in 14 C yr BP) .
σ det = The measurement uncertainty associated with A det (in 14 C yr).
A cc (t) = The 14 C activity (in 14 C yr BP) predicted by the calibration curve for a discrete age t.
σ cc (t) = The uncertainty (in 14 C yr) associated with A cc (t).

R(t) =
The predicted 14 C activity depletion (in 14 C yr) of A det relative to the calibration curve at discrete age t, due to a local reservoir effect (Stuiver et al. 1986). R(t) can be substituted with ΔR(t) in the case of a marine calibration curve.
Activity depletion due to R(t) is considered here by incorporating it into the calibration curve 14 C activity. This approach to handling R(t) allows, if desired, for temporally dynamic R(t) to be correctly incorporated (Waelbroeck et al. 2019). The calibration curve is adjusted as follows, for each discrete calendar age t: Uncertainties pertaining to calibration curve 14 C activity and the 14 C reservoir effect (σ cc (t) and σ R (t)) are both Gaussian, so they can be easily propagated into one term, for each discrete calendar age t: σ ccR (t )= √ (σ cc 2 (t)+σ R 2 (t)) Before proceeding, all of the above 14 C-related parameters are first converted into F 14 C space to facilitate more accurate calculations that take isotope mass balance into account, which is especially relevant in the case of wide range of 14 C activity (Erlenkeuser 1980;Bronk Ramsey 2008;Keigwin and Guilderson 2009), such as is the case with bioturbated sediment archives.
A sequence of probabilities can describe the closeness of a sequence of 14 C activities predicted for all discrete ages t (represented as T) available within the calibration curve (i.e. A ccR (T)), to a single 14 C activity predicted by the calibration curve for a discrete age t (i.e. A ccR (t)). This closeness, which includes a quantification of calibration curve and reservoir effect uncertainties, can be evaluated using a normal distribution for each instance of t, summing through all n values available in T to give the total relative 14 C probability for each t: )) (Eq. 6)

The prior calibration process
The prior calibration process involves moving the p prior distribution along a sliding window of calendar ages and each time computing the the hypothetical laboratory mean 14 C activity determination (h det ) that would result from each p prior placed at a sliding window starting at each t: Subsequently, it is possible to evaluate the single probability value of each h det (t) as a function of its closeness to the normal distribution of the sample's observed laboratory determination A det ± σ det : For each sliding window placed at each t, a vector of calibrated age probabilities is calculated, corresponding to each discrete age in the sliding window: This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30. Feedback is welcome, please contact the author.
Subsequently, each p cal (t) is sorted into a large matrix, referred to here as M cal (T): The final credible calibrated probability distribution corresponding to all ages T can be calculated simply by summing all rows in M cal (T): All elements in the resulting vector p cal (T) are subsequently normalised such that they sum to 1.

Script for automated calibration (biocal)
Here, a fully documented Matlab function (biocal.m) is provided for automated calculation of the procedures outlined in this study, with full compatibility in Octave. Other programming language versions of the script (e.g. Python, Julia and R) are forthcoming and will be uploaded to the same software repository upon completion. The biocal script takes full advantage of computer memory to carry out calculations using vectorised programming, thus resulting in a time-optimised routine. In the calibration protocol described here, it is assumed that it is possible to calculate P prior sliding windows along the the entire history covered by the calibration curve. However, as it would be computationally prohibitive to calibrate for the entire history of the calibration curve, biocal restricts its P prior sliding window calculations to an interval of the calibration curve covering a 3σ distance in each direction from the laboratory 14 C determination, with added padding to accommodate a long tail of P prior sitting at 3σ sigma distance. In future, when computer memory and This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30. Feedback is welcome, please contact the author.
processor power increases by another order of magnitude, it will be possible to compute sliding windows across the entire calibration curve, assuming that would ever be deemed necessary. The calculation time and memory usage is inversely related to the SAR, BD, 14 C measurement error, calibration curve uncertainty and reservoir effect uncertainty. Testing using Matlab 2020a on a Linux system with an Intel i7-9700 CPU resulted in the following times and memory usage: a Younger Dryas aged sample with SAR of 4 cm ka -1 and BD of 10 cm required 1.7 s calculation time and 2GB memory; the same sample, but with a SAR of 20 cm/ka -1 , required 0.2 s to calculate and used 100 MB of memory.

Evaluating calibration using sedimentological priors
Here, a test is carried out to determine if a the calibration protocol incorporating sediment priors results in an improved calibration process (i.e. a better estimation of the true age distribution of the measured sample) for a number of SAR scenarios using a globally representative BD of 10 cm.
First, the established understanding of bioturbation (Berger and Heath 1968) is used to calculate the associated annualised age distribution that would be expected for a discrete-depth, 1 cm sediment sample ( Fig. 1). In all cases, the mean value of the age distribution is set at 12 ka, and it is assumed that the oldest 10% of the foraminifera are broken, not picked and, therefore, not included in the distribution. This age distribution represents the ground-truth age distribution of our virtual sample, the target age distribution which would be the ideal calibrated age result. Subsequently, we can carry out a 'virtual AMS analysis' upon the ground-truth distribution by using the IntCal20 (Reimer et al. 2020) calibration curve to determine the mean 14 C activity that could be expected, in a bestcase scenario, to result from the aforementioned age distribution. For simplicity's sake, no reservoir effect is included in this demonstration, and it is assumed that the mean 14 C activity reported by IntCal20 perfectly represents the 14 C activity recorded by the sediment archive, with linear interpolation applied where necessary to achieve annual resolution.
Assuming an appropriate AMS machine error of ±80 14 C yr, the mean 14 C activity can then be calibrated in two ways, which can subsequently be compared: (1) using IntCal20 and IntCal20, supplemented by the SAR and BD priors associated with each scenario, to carry out the new calibration protocol outlined in this study.
As could be expected, the calibration protocol using sedimentological priors outperforms the standard calibration procedure in estimating the ground-truth age distribution, as shown in Fig. 1 for a number of SAR scenarios ranging between 4 and 20 cm/ka, with BD of 10 cm and constant temporal species abundance. In such use case scenarios, using the calibration protocol with sedimentological priors demonstrably leads to a more accurate calibrated age distribution, which would be ideal for improving age-depth modelling routines.
In Fig. 2, the SAR scenarios from Fig. 1 are repeated in the case of a much older ground-truth scenario (mean age of 32 ka), whereby Gaussian uncertainties associated with both the sample 14 C activity (±300 14 C yr assumed here) and the 14 C calibration curve are markedly increased. In Figs.
2e,f,g,h,i it can be seen that these larger uncertainties, when combined with increasing SAR, lead to the sedimentological priors becoming overwhelmed by the Gaussian 14 C uncertainties and, consequently, the calibrated age distribution determined by the procedure starts to approach a normal distribution. In these use case scenarios, the calibration protocol using sedimentological priors does not necessarily offer any advantage over the traditional calibration method.

Evaluating calibration using sedimentological and abundance priors
Temporal changes in species abundance (e.g. of foraminifera) will affect the shape of the species' age distribution for a given discrete depth. Here, a sine wave with a wavelength of 2000 years is used, purely for demonstrational purposes, as a theoretical temporal abundance function (Fig. 3). In Fig. 4, the same SAR scenarios as in Fig. 1 are analysed, but this time with the application of the abundance aspect. Firstly, the aforementioned sinusoidal temporal abundance function is applied to the ground truth distribution. Subsequently, the same abundance function is used as an additional prior input when running biocal, to complement the sedimentological priors. The results in Fig. 4 demonstrate how known information about temporal changes in species abundance can be used to produce better informed calibrated age estimations for bioturbated sediment archives. In Fig. 5, the 2000 year wavelength abundance function is also applied to in the case of an older ground-truth distribution (i.e. that of Fig. 3), demonstrating that abundance priors can also be used as a tool to better constrain 14 C analysis of older samples with greater uncertainty. This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30. Feedback is welcome, please contact the author.

Advice for determining prior values
In order to carry out the calibration protocol detailed here, prior values for SAR, BD, fraction broken foraminifera, temporal species abundance and temporal reservoir effect are required. A first order estimate for the sediment accumulation rate can be ascertained by examining the general relationship between age-depth determinations (including 14 C-derived age estimates based on existing calibration methods without sedimentological priors). It is possible to use an approximate prior for bioturbation depth using an estimate based on globally representative values (generally between 8 and 12 cm) (Trauth et al. 1997;Boudreau 1998). One could also directly estimate for the sediment archive itself based on 14 C investigations of the core top (Peng et al. 1979;Trauth et al. 1997;Henderiks et al. 2002), or by using 14 C measurements on single foraminifera (Lougheed et al. 2018) or, more accessibly, by measuring 14 C on a number of samples with low numbers of foraminifera and using a statistical analysis of the sample variation to infer downcore bioturbation depth ).
The fraction of unpicked, fragmented microfossils can be estimated by simply investigating the sample material (Le and Shackleton 1992). There is a risk, however, that the very oldest microfossils of the original population are completely dissolved and are therefore no longer present in the sample material as broken material (Ruddiman and Heezen 1967), which could affect assumptions regarding the p prior age distribution. In any case, one can take into account the susceptibility of a particular species to breakage (Boltovskoy 1991;Boltovskoy and Totah 1992) in combination with knowledge of bottom water chemistry (Ruddiman and Heezen 1967;Parker and Berger 1971), as well as the average residence time in the bioturbation zone, itself a function of SAR and BD .
Additional challenges are associated with determining temporal changes in species abundance, seeing as the abundance record sourced from the depth domain (i.e. the downcore, discrete-depth record) is itself modified by bioturbation (Lougheed 2020), and therefore does not reflect the species using such an approach remains a challenging task. Temporal reconstructions of abundance represent an inherent difficulty for the interpretation not just of 14 C chronological data, but downcore, multi-specimen microfossil records in general (Bard 2001;Löwemark and Grootes 2004;Löwemark et al. 2008;Lougheed 2020). A suitable approach could involve applying multiple plausible abundance scenarios when calibrating 14 C dates using the calibration protocol outlined here, and examining if the spread of calibrated age outcomes significantly affects the geochronological interpretation.

Conclusion
Current 14 C calibration routines for sediment archives do not incorporate information about sedimentological processes such as SAR and BD, meaning that current 14 C-based geochronologies systematically underestimate the total age range of a multi-specimen sample, and potentially also contain age-depth artefacts. By taking account of sedimentological processes in addition to 14 C uncertainties, a more credible calibrated age distribution can be ascertained using the protocol outlined here. It should be noted that this new calibration protocol offers most improvement in the case of lower SAR typical of deep-sea sediment archives. Determining complex, time-dependent variables such as species abundance is a complex task, and end users are encouraged to experiment with multiple, plausible abundance scenarios determine how they may affect the geochronological interpretation of their data. Similarly, it is possible to experiment with a range of plausible SAR and BD values. Gaining insight into uncertainties through this type of experimentation is important, because SAR itself can influence the age distribution (and hence 14 C activity distribution) of a sample, but in order to determine SAR accurately one needs to know the age of the sediment. This This EarthArxiv non-peer reviewed preprint was submitted to the journal Radiocarbon on 2021-06-30. Feedback is welcome, please contact the author.
Wheatcroft RA. 1992. Experimental tests for particle size-dependent bioturbation in the deep ocean.