the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 24 Jun 2026
CIR Process Age Inference Algorithm v1.0: scalable and consistent sedimentation rate modeling for ocean sediment cores via the Cox-Ingersoll-Ross process
Abstract. Deep-sea sediment cores provide a key archive of past climate, with geochemical measurements of microfossils providing both environmental information and age constraints; however, constructing continuous age-depth models is challenging due to sparse and uncertain direct dating (e.g., radiocarbon up to 50 kyr BP) and the need to integrate more densely sampled but indirect proxies such as benthic δ18O, which requires pattern alignment to reference stacks while maintaining physically plausible sediment rates. Bayesian frameworks have therefore become standard, with approaches such as BACON (Blaauw and Christen, 2011) that models inverse sedimentation rates via an autoregressive gamma process and BIGMACS (Lee et al., 2023) that integrates both direct and indirect constraints through probabilistic alignment and empirically informed priors. Despite their practical utilities, these method exhibit key limitations: BACON relies heavily on user-specified hyperparameters that are not statistically inferred from the data, while BIGMACS could employ a sedimentation-rate model defined on uneven depth grids, potentially leading to inconsistent smoothness and sensitivity to proxy resolution. Here we propose the Cox-Ingersoll-Ross (CIR) Process Age Inference Algorithm to tackle the aforementioned limitations. This multi-layer Bayesian hierarchical framework employs the CIR process as a prior on the inverse sedimentation rates to guarantee a consistent smoothness over depths to address the drawback of BIGMACS, and it allows estimation of the CIR process model parameters via the Expectation-maximization (EM) algorithm. To validate our framework, we first estimated the model parameters from a carefully curated dataset of 79 radiocarbon records, and then applied the algorithm to four other radiocarbon-dated benthic δ18O sediment core records to compare the performance to BIGMACS. The resulting age models not only show greater consistency and robustness but also preserve smoothness of posterior sedimentation rates over depths and successfully avoid alignment artifacts.
- Preprint
(3178 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
- RC1: 'Comment on egusphere-2026-2299', Anonymous Referee #1, 08 Sep 2026
-
RC2: 'Comment on egusphere-2026-2299', Anonymous Referee #2, 09 Sep 2026
This article presents an advancement of the author(s) previous age-depth modelling algorithm "BIGMACS" (Lee et al 2023). The new work is a significant advance, it uses an alternative to the AR(1) gamma process used by Bacon and BIGMACS - the Cox-Ingersoll-Ross (CIR) process - which allows for consistent smoothness along a modelled record - in contract to BIGMACS, in which smoothness varied with the local density of data points.
Like BIGMACS, the new approach can combine both direct (e.g. 14C) and indirect (e.g. matching to a reference stack) age control data, but in a more consistent way. The manuscript focusses on modelling sediment cores with radiocarbon dates, supplemented by matching to a d18O stack.
I was able to get the Python code running, but it is not really packaged sufficiently well for general use. The Github repo referenced from Zenodo was not reachable. In the Zenodo snapshot there was no licensing info and limited documentation.
I have tested the referenced Python code, and it seems to work well for a set of difficult to fit cores. As implied in the introduction sections, the CIR based algorithm is much less sensitive to user chosen meta parameters than is Bacon, particularly the resolution of the model (thickness of modelled sections in Bacon, number of inducing depth in CIR).
I have personal experience using Bacon, but not BIGMACS. In the manuscript, the performance of the new approach is compared only to the author(s) previous age-modelling algorithm "BIGMACS", perhaps for diplomatic reasons. However, I think the most widely used age-depth modelling software for radiocarbon is still Bacon. If they wish to convince people to adopt their new method it would help to have some direct comparisons.
The tunable parameters of the algorithm have been calibrated against a large set of sediment cores - but this set was restricted to marine sediment cores with relatively high sedimentation rates (> 8 cm / kyr). The four example cores also have relatively high and similar sedimentation rates. So it is unclear how well the default parametrisation would work for low sedimentation rate environments or terrestrial cores.
One further concern I have is that d18O values are supplied to the algorithm with no uncertainties, this might put too much weight on the d18O matching given that d18O values can be very uncertain.
I consider this work suitable for publication in GMD as long as the issues with the documentation of the code are corrected. I would recommend adding direct comparisons to Bacon for the work to have more impact in the field.
Detailed comments on the text.
Figures and tables
I found the figure and table legends too brief.
For figs 2,3,4,5 it is difficult to see the differences between CIR and BIGMACS - could the difference be plotted directly?
Abstract
I like the terminology "direct-" and "indirect-proxies" for age to refer to e.g. radiocarbon and d18O values that will be compared to a reference stack
line 65
"Query sedimentation records may then be standardised and integrated into the framework for age estimation."
I guess this means "new" sedimentation records for which an age-depth model is wanted?
line 147
Although technically it is an hierarchical model, I think using this term might be misleading. The description of the model is a closer to a latent process model: a latent stochastic process, deterministically transformed, then observed via a likelihood. It is not hierarchical in the sense of multilevel regression, with regularisation / shrinkage / partial pooling of the parameters via hierarchical levels.
193
"is a non-informative improper prior on the bias random variable B" - should this be biased random variable B?
271-274
Is the ability to estimate these shifts part of the existing code?
284-285
This was for the training data, but can the model handle variable reservoir ages at the inference stage?
323
typo "inducing"
326-327
"Since the record-specifcc standardisation parameter is defined as the median of the average sedimentation rates, it is not surprising that the estimated values of β are close to the corresponding α."
is this because the mean of a gamma is shape / rate and the mean has been set to ~1 by the rescaling, so alpha must be approximately beta by construction?
328-331
"Although V(0.5) is nearly identical to – albeit marginally greater than – V(2.0) for certain interval lengths, the corresponding estimate of ρ is close to zero. This would imply that sedimentation rates exhibit virtually no autocorrelation with respect to depth, a conclusion which contradicts established scientific understanding."
By this you mean you can exclude the apparent good fit at alpha = 0.5 because the fitted rho is unrealistic - why is a good fit achieved in this parameter space?
358-359
I think you mean the lower part of the core (top of the figure, but the deepest part of the core)
361-367
It is not really possible to see the differences between CIR and BIGMACS from these figures. Can you plot the difference directly?
388-394
It would be helpful to have a summarising figure showing some metric of "fit" as a function of grid resolution.
Citation: https://doi.org/10.5194/egusphere-2026-2299-RC2 -
EC1: 'Comment on egusphere-2026-2299', Heather Kim, 15 Sep 2026
Dear Authors,
Please see the following review from the third reviewer and address the comments by the deadline provided. Best regards, Heather Kim
----
The paper discusses an alternative chronology building inference process in paleo-studies. In a similar approach to other available methods, but using a CIS AR Gamma process, the authors use an informal step-wise hybrid inference process to fit the data. The paper is interesting and the examples are illustrative. However, I have major concerns regarding the overall use of the methodology, the authors conceptual understanding of the problem and the paper presentation.
1. Note that this is a journal article, not some lecture notes. I do not see the need to explain 1) Gaussian Processes, 2) mid point rule, 3) HMC and 4) SVGD. These are standard tools that only need to be mentioned and cited. A supplementary material could be used for a more pedagogical explanation.
2. The model in (10) is a piece-wise linear model, adding the mid point for additional smoothness, and the accumulations rates are modeled with a (CIS) Gamma process. This is as in Blaauw and Christen (2011), but changing the gamma process to the CIS, and the additional knot for integration. The CIS process is governed by two parameters, β and ρ, and these need to be set a priori (see l. 198). An interactive scheme is used in which the chronology is sampled from the posterior and then these parameters are updated. This hybrid scheme runs without justification, with no clear objective funtion in sight, nor a posterior distribution to sample from. What are the authors doing? This is surely not Bayesian; but indeed, “our algorithm returns the age model in the form of samples, not their single estimates only” (l. 203), but samples of what?
3. Why not run an MCMC (HMC) on the full posterior, as in Blaauw and Christen (2011)?
4. Was it worth the great added complication of the CIS model? Using the OU −Γ process is quite simpler, as in Blaauw and Christen (2011). Did the authors make a comparison?
5. The authors are worried about prior specification, but point “estimates” of hyperparameters, in some ad hoc procedure, do not assure us that priors are now estimated from data, and therefore are “better”. Input of prior knowledge of accumulation rates is critical here, given the few radiocarbon, and others, calendar markers along cores. That knowledge is simply not in data, and thus an informative, and committed, prior needs to be used to fill in the gaps.
Citation: https://doi.org/10.5194/egusphere-2026-2299-EC1
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 183 | 85 | 20 | 288 | 8 | 8 |
- HTML: 183
- PDF: 85
- XML: 20
- Total: 288
- BibTeX: 8
- EndNote: 8
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
In this article the authors present a new age modeling algorithm that uses a multi-layer Bayesian framework to infer sedimentation rates (and then age-depth models) in marine cores using δ18O and 14C data. The key advances over past work is more realistic smoothness and process-based fitting of sedimentation rates which was a limitation in BIGMACS, the previous approach to this work.
The manuscript is well written and justified, and the work is important. I have several suggestions to revision to the manuscript and codebases, but believe the paper should be published following these revisions.
I have two major suggestions for the manuscript itself.
Beyond the manuscript, I also reviewed the Matlab and python code shared as zenodo repositories. Given that a primary goal of the study is to enable others to use the methodology presented here to create age models, the functionality and useability of the code is in many ways as important as the manuscript itself. First, I tested reproducing an age model for MD95-2042 using both the Matlab and python approaches, and it worked well. The median ages were very similar using different random seeds and all of the posterior age paths are strictly monotone in depth, and the credible intervals brackets the medians. This all supports the science presented in the manuscript.
I was able to get the Matlab code to run on without statistics or parallel toolboxes in Matlab 2026a by only replacing a single function (normrnd), and could replace that with a very simple (4-line) function. The authors might consider making a similar change to support users who don't have or don't want to install those toolboxes. After making this change, the δ¹⁸O scale/shift and sedimentation rate all agree across this Matlab implementation, the Matlab results that came in zenodo, and python across all four shared cores. This all suggests that the code and results are robust, which is great.
However, a couple minor issues in the code prevented me from running them out of the box:
Both of these are simple fixes that need to be done so that people can replicate your results without fixing code.
A couple other thoughts about the code:
The other major issue is the documentation. There is no readme, license, requirements or environment files. I assume these all come from a github repository (one is mentioned on zenodo), but it is dead link (I suspect it's private). Having this code in a well described public git repository is a necessary improvement to help users find and use the code, as well as report issues and make contributions. I believe that the authors plan to add this to github, as it's referenced in a few places, and I encourage them to do this sooner than later.
There are other more minor issues with the documentation: