Bayesian age models and stacks: Combining age inferences from radiocarbon and benthic &delta;<sup>18</sup>O stratigraphic alignment

Lee, Taehee; Rand, Devin; Lisiecki, Lorraine E.; Gebbie, Geoffrey; Lawrence, Charles E.

doi:https://doi.org/10.5194/egusphere-2022-734

Preprints

https://doi.org/10.5194/egusphere-2022-734

Preprints

07 Sep 2022

| 07 Sep 2022

Bayesian age models and stacks: Combining age inferences from radiocarbon and benthic δ¹⁸O stratigraphic alignment

Taehee Lee, Devin Rand, Lorraine E. Lisiecki, Geoffrey Gebbie, and Charles E. Lawrence

Abstract. Previously developed software packages that generate probabilistic age models for ocean sediment cores are designed to use either age proxies (e.g., radiocarbon or tephra layers) or stratigraphic alignment (e.g., of benthic δ¹⁸O) and cannot combine age inferences from both techniques. Furthermore, many radiocarbon dating packages are not specifically designed for marine sediment cores and default settings may not accurately reflect the probability of sedimentation rate variability in the deep ocean, requiring subjective tuning of parameter settings. Here we present a new technique for generating Bayesian age models and stacks using ocean sediment core radiocarbon and benthic δ¹⁸O data, implemented in a software package named BIGMACS (Bayesian Inference Gaussian Process regression and Multiproxy Alignment of Continuous Signals). BIGMACS constructs multiproxy age models by combining age inferences from both radiocarbon ages and benthic δ¹⁸O stratigraphic alignment and constrains sedimentation rates using an empirically derived prior model based on 37 ¹⁴C-dated ocean sediment cores (Lin et al., 2014). BIGMACS also constructs continuous benthic δ¹⁸O stacks via a Gaussian process regression, which requires a smaller number of cores than previous stacking methods. This feature allows users to construct stacks for a region that shares a homogeneous deep water δ¹⁸O signal, while leveraging radiocarbon dates across multiple cores. Thus, BIGMACS efficiently generates local or regional stacks with smaller uncertainties in both age and δ¹⁸O than previously available techniques. We present two example regional benthic δ¹⁸O stacks and demonstrate that the multiproxy age models produced by BIGMACS are more precise than their single proxy counterparts.

Received: 11 Aug 2022 – Discussion started: 07 Sep 2022

Download & links

Preprint (PDF, 2006 KB)

Supplement (1277 KB)

Download & links

Taehee Lee et al.

Status: closed

RC1:
'Comment on egusphere-2022-734', Anonymous Referee #1, 06 Oct 2022

This paper presents an interesting evolution of the the previous d18O alignment and stacking functions. The main advancements are the the use of an empirically derived sedimentation rate prior and incorporation additional age information including radiocarbon dates and tephras, tie points, etc. While I think this will be a fine contribution to the field, I think the manuscript needs to be further clarified. The bulk of the background focuses on d18O alignment. The radiocarbon descriptions can be somewhat abbreviated and overly simplistic.

The main issue I have is that the information from which the prior was obtained is completely absent. The reader knows nothing about the sediment cores, their locations, age ranges, and depositional environments. An empirically derived prior that replaces tunable parameters is only an advancement if it is appropriate to the readers sediment core. Statements around Line 220 seem to indicate that the prior is a poor match for the data, and the prior, not the radiocarbon dates are the greatest influence on the age model in the radiocarbon-dated interval. Hence, I recommend the discussion of the prior be greatly enhanced.

I suggest the authors take a close look at spellings of acronyms and names. I noticed several different (and incorrect) spellings of Obrochta, as well as some acronyms that were transposed. I've noted the former in the below line-by-line comments.

Finally, the manuscript needs to state the system requirements for running the software. I note that it uses a parallel for loop, which is not included in the standard matlab distribution. So users without the Parallel Computing Toolbox cannot use this. Also of course the sampling of PDFs will require the stats toolbox. In my experience, most people have the stats toolbox but fewer can run a parfor loop.

***

Line by line comments

Lines 9 - 11: "...designed to use either age proxies (e.g., radiocarbon or tephra layers) or stratigraphic alignment (e.g., of benthic δ18O) and cannot combine age inferences from both techniques."

This is a bit misleading because other Bayesian models can indeed use oxygen isotope information -- just not in the way that is done in the present paper. A manually identified oxygen isotope tie point to a reference series can be used together with radiocarbon or tephra, so I suggest this first sentence be reworded.

Lines 28 - 30: "However, this method is restricted to the last 50 ka BP, suffers from variable surface reservoir ages ..., and is often low resolution causing the age model to be highly dependent on assumptions regarding sediment accumulation rate variability."

This is a somewhat outdated viewpoint. An increasing number of labs have installed AMSs, advances in automation have reduced preparation time and cost, and the ability to reliably measure radiocarbon on trace amounts of samples is making it increasingly possible to perform high resolution radiocarbon dating. Also, technically, the current marine calibration curve now extends to 55 ka, not 50 (though I agree that radiocarbon dating at those extreme ages is problematic, and resolution does not match oxygen isotope data). One study that comes to mind is:

Ishiwa, T. et al. Temporal variation in radiocarbon pathways caused by sea-level and tidal changes in the Bonaparte Gulf, northwestern Australia. Quaternary Science Reviews 266, 107079 (2021).

Lines 36 - 37: "Software packages exist to produce probabilistic age models using radiocarbon ages (Blaauw & Christen 2011; Lougheed & Obrachta, 2019), but none of these probabilistically combine age inferences from both dating techniques."

This description should be improved to clarify that the authors are referring to software packages that automatically find the optimal alignment to a reference series. The current descriptions reads as if none of the radiocarbon-centric models can use oxygen isotope tie points.

Line 46: Beyond 55 ka

Lines 48 - 49: "Sedimentation rates are realistically constrained with an empirically derived prior model rather than subjective parameter settings."

I wonder how appropriate this prior is for the possible range of sediment cores that users will inevitably throw at your model? Is it possible to specify your own priors?

Lines 65 - 66: "Radiocarbon ages must be calibrated from 14C years to calendar years with a calibration curve that accounts for changes in past atmospheric 14C production rates (Reimer et al., 2020; Heaton et al., 2020)"

This is a very simplified statement. Changes in the carbon cycle is also taken into account, and quite a lot of work has gone into better understanding changes in marine reservoir age for Marine20.

Lines 67 - 68: "The uncertainty of the calibrated age is a combination of the calibration curve uncertainty, the radiocarbon measurement uncertainty, and the marine reservoir age uncertainty."

To this list should be added 1) local reservoir age offset from the global mean, Delta R, which also has its own uncertainty and is 2) temporally variable (e.g., older water at downstream upwelling sites following AMOC slowdown, etc.)

Line 76: "LGS" should be "LSG". Also perhaps Heaton 2020 should also be cited here because Marine20 includes the BICYCLE LGS-OGCM. (see above comment starting Line 65 -- Both IntCal20 and particularly Marine20 are much more sophisticated (complicated?) than just correcting for changes in production rate.

Lines 106 and 108: "trial and error"

As with the other models described, Undatable also comes with suggestions regarding parameter selection. Both of these sentences would probably be better without "trial and error". If it took hours to converge, then that would be "trial and error", but since it takes seconds, It's more like adjusting music volume to one's desired level through instantly received feedback, which is not "trial and error". I'd suggest rewriting as:

"Its quick runtime encourages parameter tuning, based on the authors' recommendations"

and

"These parameters have large effects on the resulting age model requiring the user to decide on the most appropriate values rather than using a prior model of sedimentation rate variability."

And I also suggest that the tunable parameters in the other models be similarly discussed. As it is, this description reads as is undatable is the only one with parameters that can be tuned. This somewhat undersells what the authors are presenting here: a model without tunable parameters.

Line 112: "which often correlates with salinity."

yes it might loosely correspond to salinity but it's really surface evaporation - precipitation prior to deepwater formation (since I assume the author's mainly considering benthic oxygen isotopes.

Line 115: "The most conservative technique for aligning records to a target is to assume that large, easily identifiable features in the signals, such as glacial terminations, occurred simultaneously, create tie points between these features, and linearly interpolate between the tie points"

There absolutely is a lag between "upstream" sites in the North Atlantic and "downstream" sites since it can take on the order of 1000 years or more for the signal to propagate with the flow of deepwater.

Line 135: I suggest this be better presented with the information starting Line 115.

Line 165: "termed the likelihoods" remove "the"

Line 171 - 175 The sedimentation model is called a prior distribution which is in turn called a transition model. Perhaps this can be made more clear.

Line 180: "confidence" should be "credible"?

Line 190: When are the locations of these 37 cores going to be disclosed?

Line 191: "However, where the previous study interpolated sedimentation rates every 1 kyr, we interpolate by 1 cm"

What is the range of sedimentation rates in the 37 cores? Is 1 cm sampling typically equivalent to a 1 ky sampling, or is the interpolation interval vastly different than that used by Lin et al?

Line 201 - 202: "Expansion specifies a below average sedimentation rate and refers to a stretching of the local portion of the record."

This is a bit confusing as stated and doesn't become clear until the next sentence where the authors stake that "contraction ... requires squeezing" Maybe rewrite as:

"Expansion refers to a below average sedimentation necessitation stretching the local portion of the record"

Lines 203 - 204: "If the local sedimentation rate is within 8% of the core’s average, the state is classified as steady."

How was 8% selected? Please further clarify as is done on lines 209 - 210 regarding the 15 cm interval.

Line 220: "improves agreement between the core age models and the radiocarbon observations"

I don't understand this sentence. The age model should be based on the radiocarbon observation in the radiocarbon-dated intervals. Does this indicate that the prior is often vastly different from the data, and without changing the alpha and beta parameters relative to the previous Bayesian models, the age model obtained by BIGMACS is inconsistent with the radiocarbon dates?

Lines 234 - 237: "Specifying the model as a uniform distribution will force the age model to pass through the given uncertainty range and should be used when the user is confident about the age information. Specifying a Gaussian distribution will allow the age model to pass farther from the

additional age constraint."

This seems backwards to me. If I specify a tephra age as a gaussian distribution with some mean and standard deviation, the highest probability is at the mean, so the model should pass closest to the mean. But if I specify a uniform distribution, the model has an equal probability of passing anywhere. So wouldn't the user want to specify a gaussian when there is good confidence in the age constraint? Perhaps I'm not following what the authors mean to say. Is it that when there is confidence in the *other* age data, with less confidence in the specified tephra/tie point, that the authors are suggesting to use a uniform distribution? I think this statements needs to be clarified.

Lines 301 - 302: "these cores contain a relatively large number of δ18O outliers (Figure 1)."

Not an appropriate text location to reference fig 1. Please add lat and lon to fig 1.

Table 1: confirm the longitudes

Figure 2: a color bar for the panels A and B would be helpful.

Line 347: "... crosses in Figures 4A and 5A ..."

Figure 3 has yet to be mentioned. Confirm figure numberings. I think this should be Fig 3A and 4A. Generally Figures are numbered in the order they are mentioned in the text.

Line 356: "Figure 6 compares the DNEA and ITWA stacks"

Change to Fig 5.

Line 370: "The Gaussian process regression also creates smoother stacks than previous binning methods"

It would be very useful to the reader to see a comparison of the previous stacking methods. It would also be very helpful to add a figure showing each sediment core's d18O record plotted in a separate panel above the BIGMACS stack. This will let the reader better visualize the the smoothing due to the increased autocorrelation. This would also support the assertion on Line 385 of homogenous signals.

Line 449: "6.1.1 Radiocarbon and multiproxy age models"

Missing from the discussion of applicability is, of course, if the goal is to compare phasing between d18O records, then the multiproxy age model cannot be used and only 14C, tephras, etc. can be used.

Lines 456 - 458: "Because BIGMACS applies a prior model based on observed sedimentation rate variability (Lin et al., 2014), the age uncertainty between 14C observations returned by BIGMACS is physically realistic and less subjective than using tuned parameters in other software packages."

At this point, we still know nothing about the cores from which this prior was obtained. Where are they located? What are their water depths? What are their age ranges? Do they span glacial/interglacial terminations? While this methods does not require parameter selection, it is assuming that the prior is reasonable for the *user's* sediment cores. This is an extremely important point, and I think the authors should spend some time to demonstrate to the reader that the prior is actually appropriate. In short, I'd like to have it explained to me very clearly why the prior assumed here is both appropriate and better than selecting parameters. The statement I mentioned earlier on Line 220 gives the impression that the prior is overly informative and inconsistent with the data.

Line 470: "widely space"

"spaced"

Lines 478 - 479: "an assessment of a core’s absolute age uncertainty should incorporate both the absolute age uncertainty of the target/stack and alignment uncertainty."

I would suggest adding an optional age error column for the stacking target, then fold that error into the alignment uncertainty. You could output both age uncertainty obtained from that of the alignment target, in addition to the alignment uncertainty already returned. The could be added to get a total uncertainty.

Line 544: "Example multiproxy regional stacks"

The age models are "multiproxy" but the stacks are not.

Line 547: "standard deviations include the effects of spatial variability, age uncertainty ..."

I really think that there should be an easier way for users to include the age uncertainty in the alignment target.

Author contributions: It appears that the first two listed authors contributed equally. As such they should be listed as "contributed equally" somewhere around where the corresponding authors are noted. If the other authors only contributed funding for this study, then technically they should not be authors and should be acknowledged.

***

I didn't do a full code review but I do have some suggestions as the authors suggest that BIGMACS is resource intensive and slow. There are several things that I see that could be optimized. While I feel that the time and memory savings on the things I am point out will be minimal, it makes me wonder if there are similar inefficiencies in the most critical parts of code.

getInitialTarget.m

Line 66, load calibration curve

Why not just load only the curves that are needed? There are much more efficient ways to read in the data. The fastest is just remove the headerline of each calcurve and use simply load(path). Small things like this, if they occur throughout the code base, can add up to a savings in runtime. Also note that "path" is a command to Get/set search the path. I'd suggest changing the variable name to "Path" or "pth".

Why do all this:

tic

path = 'Defaults/Calibration_Curves/IntCal20.txt';

fileID = fopen(path);

CAL = textscan(fileID,'%s %s %s %s %s');

fclose(fileID);

cal_curve{1} = zeros(length(CAL{1})-1,5);

for k = 1:5

cal_curve{1}(:,k) = str2double(CAL{k}(2:end));

end

toc

Elapsed time is 0.390675 seconds.

when you can do simply this, which is simpler and an order of magnitude faster. Are there similar chunks of inefficient code that are resulting in slow runtime?

tic

Path = 'Defaults/Calibration_Curves/IntCal20.txt';

fileID = fopen(Path);

CAL = textscan(fileID,'%d %d %d %d %d','headerlines',1);

toc

Elapsed time is 0.018334 seconds.

getData.m

If you can figure out what the final size of e.g., "d18O_depth" will be, you can preallocate a matrix for better memory management and faster runtime.

initializeAlignment.m getAlignment.m

This function uses a parfor loop, which requires the parallel computing toolbox that not everyone will have. It will also take time to start up the parallel pool if not already running. Could check for the existence of the toolbox and if it's not installed, use for instead. If there is not a significant improvement in speed, giving the time to start up the pool, it might be better to use just a for loop.

Citation: https://doi.org/10.5194/egusphere-2022-734-RC1
- AC1: 'Reply on RC1', Devin Rand, 08 Feb 2023
  
  The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2022/egusphere-2022-734/egusphere-2022-734-AC1-supplement.pdf
  
  Citation: https://doi.org/10.5194/egusphere-2022-734-AC1
RC2:
'Comment on egusphere-2022-734', Tim Heaton, 28 Dec 2022

Thanks to the authors for their work and letting me contriibute to the discussion. Please see attached file for my thoughts,

Tim Heaton

Citation: https://doi.org/10.5194/egusphere-2022-734-RC2
- AC2: 'Reply on RC2', Devin Rand, 08 Feb 2023
  
  The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2022/egusphere-2022-734/egusphere-2022-734-AC2-supplement.pdf
  
  Citation: https://doi.org/10.5194/egusphere-2022-734-AC2
EC1:
'Comment on egusphere-2022-734', Luke Skinner, 10 Jan 2023

Now that a sufficient number of reviews have been submitted, I would like to invite you provide a point-by-point response to each of the review comments. I would like to underline the importance of addressing a key issue that has been raised in the reviews, namely: the ‘empirically constrained’ sedimentation rate prior that is applied in the matching algorithm. One issue is that the validity and applicability of this prior, across a range of sedimentary contexts, does not appear to have been fully assessed in a transparent manner in the manuscript – and indeed seems doubtful.

Below I add a few remarks of my own, in case these are helpful for considering what revisions you would undertake, and I invite your response to these too. I find your study of particular interest, and I hope that my comments will be seen as useful.

Title/general ethos:
In general, I think it might be useful to more clearly delineate the distinction between alignment and ‘dating’ at the core of the manuscript (even though the difference between relative and ‘absolute/numerical’ ages is indeed noted in the paper a few times). Creating a benthic d18O stack is one thing, aligning to a benthic d18O stack is another, and dating a sediment core is yet another. The only way that benthic alignment provides age constraints is if one proposes to have prior knowledge of how local/regional deep-water T and d18Osw relate to insolation, e.g. based on a hypothesis for how insolation paces ice volume, and how changes in ice volume are linked to deep water T changes and/or influence deep ocean d18Osw at a given location in the ocean. The latter sequence of hypotheses can give age constraints that are of ~millennial accuracy at best. In such a context, radiocarbon dates (even with ~centennial uncertainties in reservoir age offsets) obviously can provide a refinement. The inverse is unlikely to be true: age constraints on benthic d18O are unlikely to be precise enough, even to constrain changes in radiocarbon reservoir age offsets of order 100-1000 years. Alternatively, if the core notion of the manuscript and algorithm is the simple transferral of a radiocarbon chronology (or the pooling of radiocarbon dates) between sediment cores via a stratigraphic alignment of benthic d18O, then again it is not quite a case of ‘combining age constraints from radiocarbon and benthic d18O’. Rather, it is one of radiocarbon dating of a stratigraphic alignment/stack.

As such, my own feeling is that the manuscript might more accurately be framed in terms of ‘refining orbitally-tuned benthic d18O age models using radiocarbon constraints’, e.g. in the title and through the text. In a similar vein, it seems to me that describing the age models as ‘multi-proxy’ is a little misleading: my own expectation was initially of something like that described in line 522. I would again suggest that the process tackled in the present study be described as something like ‘radiocarbon-refined single proxy alignment’.

Line 25:
This line is not quite correct: the accuracy with which ocean sediment cores can reconstruct the *timing* of past climate events, depends on.. the.. age model. The accuracy of proxies is a separate (thorny) matter.

Line 48:
“Sedimentation rates are realistically constrained….”
As pointed out by the first Reviewer, it seems we must take this on faith, whereas there is burden of demonstration here.

Line 65:
In general, there is a need to be precise when describing radiocarbon procedures. Radiocarbon dates need to be calibrated to account for past changes in the initial radiocarbon concentration of the fossil entity’s ‘parent reservoir’ (atmosphere, surface ocean, etc.), which may change due to 14C-production changes and/or other carbon cycle processes. This crops up again on Line 72: planktonic foraminiferal radiocarbon dates must be corrected for ‘reservoir age offsets’ (relative to the atmosphere) only if using a record of past atmospheric radiocarbon concentration/activity for the calibration. In principle, a ‘marine calibration curve’ might be used instead, with different potential corrections needed as a result.

Line 79:
“…requires simulating the core’s sedimentation rate.”
I think this might be more accurately phrased as: “…requires the assumptions/models of the core’s evolving sedimentation rate between dated intervals.”

Line 90:
I think this is a but unfair to Bchron: instead of ‘resulting in extreme sedimentation rate variability’, it simply posits the full range of possibility wherever there are no prior constraints on sedimentation rates. This is arguably pretty sensible, and it represents a useful counter point to methods that assume a priori knowledge of sedimentation rates.

Line 109:
Again on the sedimentation rate prior issue: does a prior on sedimentation rate not ‘beg the question’ with regard to down-core changes in age, requiring simply a single point to be anchored in time? This seems like a very (overly) strong constraint to apply, does it not?

Line 138:
Is it worth noting perhaps that this shifts the problem of assuming ‘instant ocean mixing’ to one of a priori knowledge of past ocean hydrography and circulation?

Table 1:
note that the 14C dates for MD99-2334K are reported only by Skinner et al., G-cubed 2003 (Skinner & Elderfield 2003 does not exist, and was omitted from the references for this reason no doubt); Skinner and Shackleton 2004, and Skinner et al., Paleoc. & Paleoclim. 2021.

Figure 6: What is the reason for choosing this sediment core in particular? MD99-2334K is included in the present study, has various alternative stratigraphic age-models (aligned to the Greenland ice core event stratigraphy, and the Hulu speleothem record), as well as a reasonable 14C chronology, and a well resolved benthic d18O record. Would this not be an optimal target for testing the method? A comparison with MD95-2042 could also be made, since both also have ‘alignable’ planktic d18O records. Furthermore, these two cores were obtained using different coring devices resulting in very different ‘apparent sedimentation rates’ (due to compaction in the Kasten core and stretching in the Calypso corer), providing a useful basis for assessing the algorithm’s sedimentation rate prior.

Line 537: again, I would propose that it might be more transparent to refer to ‘radiocarbon-refined/guided d18O alignments, or similar. I wonder what the authors think.

I look forward to reading your views on these, and most importantly the reviewers’, comments.

Sincerely
Luke Skinner

Citation: https://doi.org/10.5194/egusphere-2022-734-EC1
- AC3: 'Reply on EC1', Devin Rand, 08 Feb 2023
  
  The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2022/egusphere-2022-734/egusphere-2022-734-AC3-supplement.pdf
  
  Citation: https://doi.org/10.5194/egusphere-2022-734-AC3

Status: closed

RC1:
'Comment on egusphere-2022-734', Anonymous Referee #1, 06 Oct 2022

This paper presents an interesting evolution of the the previous d18O alignment and stacking functions. The main advancements are the the use of an empirically derived sedimentation rate prior and incorporation additional age information including radiocarbon dates and tephras, tie points, etc. While I think this will be a fine contribution to the field, I think the manuscript needs to be further clarified. The bulk of the background focuses on d18O alignment. The radiocarbon descriptions can be somewhat abbreviated and overly simplistic.

The main issue I have is that the information from which the prior was obtained is completely absent. The reader knows nothing about the sediment cores, their locations, age ranges, and depositional environments. An empirically derived prior that replaces tunable parameters is only an advancement if it is appropriate to the readers sediment core. Statements around Line 220 seem to indicate that the prior is a poor match for the data, and the prior, not the radiocarbon dates are the greatest influence on the age model in the radiocarbon-dated interval. Hence, I recommend the discussion of the prior be greatly enhanced.

I suggest the authors take a close look at spellings of acronyms and names. I noticed several different (and incorrect) spellings of Obrochta, as well as some acronyms that were transposed. I've noted the former in the below line-by-line comments.

Finally, the manuscript needs to state the system requirements for running the software. I note that it uses a parallel for loop, which is not included in the standard matlab distribution. So users without the Parallel Computing Toolbox cannot use this. Also of course the sampling of PDFs will require the stats toolbox. In my experience, most people have the stats toolbox but fewer can run a parfor loop.

***

Line by line comments

Lines 9 - 11: "...designed to use either age proxies (e.g., radiocarbon or tephra layers) or stratigraphic alignment (e.g., of benthic δ18O) and cannot combine age inferences from both techniques."

This is a bit misleading because other Bayesian models can indeed use oxygen isotope information -- just not in the way that is done in the present paper. A manually identified oxygen isotope tie point to a reference series can be used together with radiocarbon or tephra, so I suggest this first sentence be reworded.

Lines 28 - 30: "However, this method is restricted to the last 50 ka BP, suffers from variable surface reservoir ages ..., and is often low resolution causing the age model to be highly dependent on assumptions regarding sediment accumulation rate variability."

This is a somewhat outdated viewpoint. An increasing number of labs have installed AMSs, advances in automation have reduced preparation time and cost, and the ability to reliably measure radiocarbon on trace amounts of samples is making it increasingly possible to perform high resolution radiocarbon dating. Also, technically, the current marine calibration curve now extends to 55 ka, not 50 (though I agree that radiocarbon dating at those extreme ages is problematic, and resolution does not match oxygen isotope data). One study that comes to mind is:

Ishiwa, T. et al. Temporal variation in radiocarbon pathways caused by sea-level and tidal changes in the Bonaparte Gulf, northwestern Australia. Quaternary Science Reviews 266, 107079 (2021).

Lines 36 - 37: "Software packages exist to produce probabilistic age models using radiocarbon ages (Blaauw & Christen 2011; Lougheed & Obrachta, 2019), but none of these probabilistically combine age inferences from both dating techniques."

This description should be improved to clarify that the authors are referring to software packages that automatically find the optimal alignment to a reference series. The current descriptions reads as if none of the radiocarbon-centric models can use oxygen isotope tie points.

Line 46: Beyond 55 ka

Lines 48 - 49: "Sedimentation rates are realistically constrained with an empirically derived prior model rather than subjective parameter settings."

I wonder how appropriate this prior is for the possible range of sediment cores that users will inevitably throw at your model? Is it possible to specify your own priors?

Lines 65 - 66: "Radiocarbon ages must be calibrated from 14C years to calendar years with a calibration curve that accounts for changes in past atmospheric 14C production rates (Reimer et al., 2020; Heaton et al., 2020)"

This is a very simplified statement. Changes in the carbon cycle is also taken into account, and quite a lot of work has gone into better understanding changes in marine reservoir age for Marine20.

Lines 67 - 68: "The uncertainty of the calibrated age is a combination of the calibration curve uncertainty, the radiocarbon measurement uncertainty, and the marine reservoir age uncertainty."

To this list should be added 1) local reservoir age offset from the global mean, Delta R, which also has its own uncertainty and is 2) temporally variable (e.g., older water at downstream upwelling sites following AMOC slowdown, etc.)

Line 76: "LGS" should be "LSG". Also perhaps Heaton 2020 should also be cited here because Marine20 includes the BICYCLE LGS-OGCM. (see above comment starting Line 65 -- Both IntCal20 and particularly Marine20 are much more sophisticated (complicated?) than just correcting for changes in production rate.

Lines 106 and 108: "trial and error"

As with the other models described, Undatable also comes with suggestions regarding parameter selection. Both of these sentences would probably be better without "trial and error". If it took hours to converge, then that would be "trial and error", but since it takes seconds, It's more like adjusting music volume to one's desired level through instantly received feedback, which is not "trial and error". I'd suggest rewriting as:

"Its quick runtime encourages parameter tuning, based on the authors' recommendations"

and

"These parameters have large effects on the resulting age model requiring the user to decide on the most appropriate values rather than using a prior model of sedimentation rate variability."

And I also suggest that the tunable parameters in the other models be similarly discussed. As it is, this description reads as is undatable is the only one with parameters that can be tuned. This somewhat undersells what the authors are presenting here: a model without tunable parameters.

Line 112: "which often correlates with salinity."

yes it might loosely correspond to salinity but it's really surface evaporation - precipitation prior to deepwater formation (since I assume the author's mainly considering benthic oxygen isotopes.

Line 115: "The most conservative technique for aligning records to a target is to assume that large, easily identifiable features in the signals, such as glacial terminations, occurred simultaneously, create tie points between these features, and linearly interpolate between the tie points"

There absolutely is a lag between "upstream" sites in the North Atlantic and "downstream" sites since it can take on the order of 1000 years or more for the signal to propagate with the flow of deepwater.

Line 135: I suggest this be better presented with the information starting Line 115.

Line 165: "termed the likelihoods" remove "the"

Line 171 - 175 The sedimentation model is called a prior distribution which is in turn called a transition model. Perhaps this can be made more clear.

Line 180: "confidence" should be "credible"?

Line 190: When are the locations of these 37 cores going to be disclosed?

Line 191: "However, where the previous study interpolated sedimentation rates every 1 kyr, we interpolate by 1 cm"

What is the range of sedimentation rates in the 37 cores? Is 1 cm sampling typically equivalent to a 1 ky sampling, or is the interpolation interval vastly different than that used by Lin et al?

Line 201 - 202: "Expansion specifies a below average sedimentation rate and refers to a stretching of the local portion of the record."

This is a bit confusing as stated and doesn't become clear until the next sentence where the authors stake that "contraction ... requires squeezing" Maybe rewrite as:

"Expansion refers to a below average sedimentation necessitation stretching the local portion of the record"

Lines 203 - 204: "If the local sedimentation rate is within 8% of the core’s average, the state is classified as steady."

How was 8% selected? Please further clarify as is done on lines 209 - 210 regarding the 15 cm interval.

Line 220: "improves agreement between the core age models and the radiocarbon observations"

I don't understand this sentence. The age model should be based on the radiocarbon observation in the radiocarbon-dated intervals. Does this indicate that the prior is often vastly different from the data, and without changing the alpha and beta parameters relative to the previous Bayesian models, the age model obtained by BIGMACS is inconsistent with the radiocarbon dates?

Lines 234 - 237: "Specifying the model as a uniform distribution will force the age model to pass through the given uncertainty range and should be used when the user is confident about the age information. Specifying a Gaussian distribution will allow the age model to pass farther from the

additional age constraint."

This seems backwards to me. If I specify a tephra age as a gaussian distribution with some mean and standard deviation, the highest probability is at the mean, so the model should pass closest to the mean. But if I specify a uniform distribution, the model has an equal probability of passing anywhere. So wouldn't the user want to specify a gaussian when there is good confidence in the age constraint? Perhaps I'm not following what the authors mean to say. Is it that when there is confidence in the *other* age data, with less confidence in the specified tephra/tie point, that the authors are suggesting to use a uniform distribution? I think this statements needs to be clarified.

Lines 301 - 302: "these cores contain a relatively large number of δ18O outliers (Figure 1)."

Not an appropriate text location to reference fig 1. Please add lat and lon to fig 1.

Table 1: confirm the longitudes

Figure 2: a color bar for the panels A and B would be helpful.

Line 347: "... crosses in Figures 4A and 5A ..."

Figure 3 has yet to be mentioned. Confirm figure numberings. I think this should be Fig 3A and 4A. Generally Figures are numbered in the order they are mentioned in the text.

Line 356: "Figure 6 compares the DNEA and ITWA stacks"

Change to Fig 5.

Line 370: "The Gaussian process regression also creates smoother stacks than previous binning methods"

It would be very useful to the reader to see a comparison of the previous stacking methods. It would also be very helpful to add a figure showing each sediment core's d18O record plotted in a separate panel above the BIGMACS stack. This will let the reader better visualize the the smoothing due to the increased autocorrelation. This would also support the assertion on Line 385 of homogenous signals.

Line 449: "6.1.1 Radiocarbon and multiproxy age models"

Missing from the discussion of applicability is, of course, if the goal is to compare phasing between d18O records, then the multiproxy age model cannot be used and only 14C, tephras, etc. can be used.

Lines 456 - 458: "Because BIGMACS applies a prior model based on observed sedimentation rate variability (Lin et al., 2014), the age uncertainty between 14C observations returned by BIGMACS is physically realistic and less subjective than using tuned parameters in other software packages."

At this point, we still know nothing about the cores from which this prior was obtained. Where are they located? What are their water depths? What are their age ranges? Do they span glacial/interglacial terminations? While this methods does not require parameter selection, it is assuming that the prior is reasonable for the *user's* sediment cores. This is an extremely important point, and I think the authors should spend some time to demonstrate to the reader that the prior is actually appropriate. In short, I'd like to have it explained to me very clearly why the prior assumed here is both appropriate and better than selecting parameters. The statement I mentioned earlier on Line 220 gives the impression that the prior is overly informative and inconsistent with the data.

Line 470: "widely space"

"spaced"

Lines 478 - 479: "an assessment of a core’s absolute age uncertainty should incorporate both the absolute age uncertainty of the target/stack and alignment uncertainty."

I would suggest adding an optional age error column for the stacking target, then fold that error into the alignment uncertainty. You could output both age uncertainty obtained from that of the alignment target, in addition to the alignment uncertainty already returned. The could be added to get a total uncertainty.

Line 544: "Example multiproxy regional stacks"

The age models are "multiproxy" but the stacks are not.

Line 547: "standard deviations include the effects of spatial variability, age uncertainty ..."

I really think that there should be an easier way for users to include the age uncertainty in the alignment target.

Author contributions: It appears that the first two listed authors contributed equally. As such they should be listed as "contributed equally" somewhere around where the corresponding authors are noted. If the other authors only contributed funding for this study, then technically they should not be authors and should be acknowledged.

***

I didn't do a full code review but I do have some suggestions as the authors suggest that BIGMACS is resource intensive and slow. There are several things that I see that could be optimized. While I feel that the time and memory savings on the things I am point out will be minimal, it makes me wonder if there are similar inefficiencies in the most critical parts of code.

getInitialTarget.m

Line 66, load calibration curve

Why not just load only the curves that are needed? There are much more efficient ways to read in the data. The fastest is just remove the headerline of each calcurve and use simply load(path). Small things like this, if they occur throughout the code base, can add up to a savings in runtime. Also note that "path" is a command to Get/set search the path. I'd suggest changing the variable name to "Path" or "pth".

Why do all this:

tic

path = 'Defaults/Calibration_Curves/IntCal20.txt';

fileID = fopen(path);

CAL = textscan(fileID,'%s %s %s %s %s');

fclose(fileID);

cal_curve{1} = zeros(length(CAL{1})-1,5);

for k = 1:5

cal_curve{1}(:,k) = str2double(CAL{k}(2:end));

end

toc

Elapsed time is 0.390675 seconds.

when you can do simply this, which is simpler and an order of magnitude faster. Are there similar chunks of inefficient code that are resulting in slow runtime?

tic

Path = 'Defaults/Calibration_Curves/IntCal20.txt';

fileID = fopen(Path);

CAL = textscan(fileID,'%d %d %d %d %d','headerlines',1);

toc

Elapsed time is 0.018334 seconds.

getData.m

If you can figure out what the final size of e.g., "d18O_depth" will be, you can preallocate a matrix for better memory management and faster runtime.

initializeAlignment.m getAlignment.m

This function uses a parfor loop, which requires the parallel computing toolbox that not everyone will have. It will also take time to start up the parallel pool if not already running. Could check for the existence of the toolbox and if it's not installed, use for instead. If there is not a significant improvement in speed, giving the time to start up the pool, it might be better to use just a for loop.

Citation: https://doi.org/10.5194/egusphere-2022-734-RC1
- AC1: 'Reply on RC1', Devin Rand, 08 Feb 2023
  
  The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2022/egusphere-2022-734/egusphere-2022-734-AC1-supplement.pdf
  
  Citation: https://doi.org/10.5194/egusphere-2022-734-AC1
RC2:
'Comment on egusphere-2022-734', Tim Heaton, 28 Dec 2022

Thanks to the authors for their work and letting me contriibute to the discussion. Please see attached file for my thoughts,

Tim Heaton

Citation: https://doi.org/10.5194/egusphere-2022-734-RC2
- AC2: 'Reply on RC2', Devin Rand, 08 Feb 2023
  
  The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2022/egusphere-2022-734/egusphere-2022-734-AC2-supplement.pdf
  
  Citation: https://doi.org/10.5194/egusphere-2022-734-AC2
EC1:
'Comment on egusphere-2022-734', Luke Skinner, 10 Jan 2023

Now that a sufficient number of reviews have been submitted, I would like to invite you provide a point-by-point response to each of the review comments. I would like to underline the importance of addressing a key issue that has been raised in the reviews, namely: the ‘empirically constrained’ sedimentation rate prior that is applied in the matching algorithm. One issue is that the validity and applicability of this prior, across a range of sedimentary contexts, does not appear to have been fully assessed in a transparent manner in the manuscript – and indeed seems doubtful.

Below I add a few remarks of my own, in case these are helpful for considering what revisions you would undertake, and I invite your response to these too. I find your study of particular interest, and I hope that my comments will be seen as useful.

Title/general ethos:
In general, I think it might be useful to more clearly delineate the distinction between alignment and ‘dating’ at the core of the manuscript (even though the difference between relative and ‘absolute/numerical’ ages is indeed noted in the paper a few times). Creating a benthic d18O stack is one thing, aligning to a benthic d18O stack is another, and dating a sediment core is yet another. The only way that benthic alignment provides age constraints is if one proposes to have prior knowledge of how local/regional deep-water T and d18Osw relate to insolation, e.g. based on a hypothesis for how insolation paces ice volume, and how changes in ice volume are linked to deep water T changes and/or influence deep ocean d18Osw at a given location in the ocean. The latter sequence of hypotheses can give age constraints that are of ~millennial accuracy at best. In such a context, radiocarbon dates (even with ~centennial uncertainties in reservoir age offsets) obviously can provide a refinement. The inverse is unlikely to be true: age constraints on benthic d18O are unlikely to be precise enough, even to constrain changes in radiocarbon reservoir age offsets of order 100-1000 years. Alternatively, if the core notion of the manuscript and algorithm is the simple transferral of a radiocarbon chronology (or the pooling of radiocarbon dates) between sediment cores via a stratigraphic alignment of benthic d18O, then again it is not quite a case of ‘combining age constraints from radiocarbon and benthic d18O’. Rather, it is one of radiocarbon dating of a stratigraphic alignment/stack.

As such, my own feeling is that the manuscript might more accurately be framed in terms of ‘refining orbitally-tuned benthic d18O age models using radiocarbon constraints’, e.g. in the title and through the text. In a similar vein, it seems to me that describing the age models as ‘multi-proxy’ is a little misleading: my own expectation was initially of something like that described in line 522. I would again suggest that the process tackled in the present study be described as something like ‘radiocarbon-refined single proxy alignment’.

Line 25:
This line is not quite correct: the accuracy with which ocean sediment cores can reconstruct the *timing* of past climate events, depends on.. the.. age model. The accuracy of proxies is a separate (thorny) matter.

Line 48:
“Sedimentation rates are realistically constrained….”
As pointed out by the first Reviewer, it seems we must take this on faith, whereas there is burden of demonstration here.

Line 65:
In general, there is a need to be precise when describing radiocarbon procedures. Radiocarbon dates need to be calibrated to account for past changes in the initial radiocarbon concentration of the fossil entity’s ‘parent reservoir’ (atmosphere, surface ocean, etc.), which may change due to 14C-production changes and/or other carbon cycle processes. This crops up again on Line 72: planktonic foraminiferal radiocarbon dates must be corrected for ‘reservoir age offsets’ (relative to the atmosphere) only if using a record of past atmospheric radiocarbon concentration/activity for the calibration. In principle, a ‘marine calibration curve’ might be used instead, with different potential corrections needed as a result.

Line 79:
“…requires simulating the core’s sedimentation rate.”
I think this might be more accurately phrased as: “…requires the assumptions/models of the core’s evolving sedimentation rate between dated intervals.”

Line 90:
I think this is a but unfair to Bchron: instead of ‘resulting in extreme sedimentation rate variability’, it simply posits the full range of possibility wherever there are no prior constraints on sedimentation rates. This is arguably pretty sensible, and it represents a useful counter point to methods that assume a priori knowledge of sedimentation rates.

Line 109:
Again on the sedimentation rate prior issue: does a prior on sedimentation rate not ‘beg the question’ with regard to down-core changes in age, requiring simply a single point to be anchored in time? This seems like a very (overly) strong constraint to apply, does it not?

Line 138:
Is it worth noting perhaps that this shifts the problem of assuming ‘instant ocean mixing’ to one of a priori knowledge of past ocean hydrography and circulation?

Table 1:
note that the 14C dates for MD99-2334K are reported only by Skinner et al., G-cubed 2003 (Skinner & Elderfield 2003 does not exist, and was omitted from the references for this reason no doubt); Skinner and Shackleton 2004, and Skinner et al., Paleoc. & Paleoclim. 2021.

Figure 6: What is the reason for choosing this sediment core in particular? MD99-2334K is included in the present study, has various alternative stratigraphic age-models (aligned to the Greenland ice core event stratigraphy, and the Hulu speleothem record), as well as a reasonable 14C chronology, and a well resolved benthic d18O record. Would this not be an optimal target for testing the method? A comparison with MD95-2042 could also be made, since both also have ‘alignable’ planktic d18O records. Furthermore, these two cores were obtained using different coring devices resulting in very different ‘apparent sedimentation rates’ (due to compaction in the Kasten core and stretching in the Calypso corer), providing a useful basis for assessing the algorithm’s sedimentation rate prior.

Line 537: again, I would propose that it might be more transparent to refer to ‘radiocarbon-refined/guided d18O alignments, or similar. I wonder what the authors think.

I look forward to reading your views on these, and most importantly the reviewers’, comments.

Sincerely
Luke Skinner

Citation: https://doi.org/10.5194/egusphere-2022-734-EC1
- AC3: 'Reply on EC1', Devin Rand, 08 Feb 2023
  
  The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2022/egusphere-2022-734/egusphere-2022-734-AC3-supplement.pdf
  
  Citation: https://doi.org/10.5194/egusphere-2022-734-AC3

Taehee Lee et al.

Supplement

https://doi.org/10.5194/egusphere-2022-734-supplement

Model code and software

BIGMACS software package Taehee Lee, Devin Rand, Lorraine E. Lisiecki, Geoffrey Gebbie, Charles E. Lawrence https://github.com/eilion/BIGMACS

Taehee Lee et al.

Viewed

Total article views: 1,147 (including HTML, PDF, and XML)

HTML	PDF	XML	Total	Supplement	BibTeX	EndNote
761	360	26	1,147	120	15	15

HTML: 761
PDF: 360
XML: 26
Total: 1,147
Supplement: 120
BibTeX: 15
EndNote: 15

Views and downloads (calculated since 07 Sep 2022)

Month	HTML	PDF	XML	Total
Sep 2022	187	82	6	275
Oct 2022	87	30	2	119
Nov 2022	52	26	0	78
Dec 2022	55	24	2	81
Jan 2023	30	19	3	52
Feb 2023	68	30	6	104
Mar 2023	35	19	1	55
Apr 2023	23	16	0	39
May 2023	26	19	0	45
Jun 2023	30	14	1	45
Jul 2023	57	26	1	84
Aug 2023	68	30	1	99
Sep 2023	43	25	3	71

Cumulative views and downloads (calculated since 07 Sep 2022)

Month	HTML	PDF	XML	Total
Sep 2022	187	82	6	275
Oct 2022	87	30	2	119
Nov 2022	52	26	0	78
Dec 2022	55	24	2	81
Jan 2023	30	19	3	52
Feb 2023	68	30	6	104
Mar 2023	35	19	1	55
Apr 2023	23	16	0	39
May 2023	26	19	0	45
Jun 2023	30	14	1	45
Jul 2023	57	26	1	84
Aug 2023	68	30	1	99
Sep 2023	43	25	3	71

Viewed (geographical distribution)

Total article views: 1,133 (including HTML, PDF, and XML) Thereof 1,133 with geography defined and 0 with unknown origin.

Country	#	Views	%

Cited

Latest update: 30 Sep 2023

Short summary

Understanding of past climate change depends, in part, on how accurately we can estimate the ages of events recorded in geologic archives. Here we present a new software package, called BIGMACS, to improve age estimates for paleoclimate data from ocean sediment cores. BIGMACS creates “multiproxy” age estimates that reduce age uncertainty by probabilistically combining information from direct age estimates, such as radiocarbon dates, and alignment of regional paleoclimate time series.


Total:	0
HTML:	0
PDF:	0
XML:	0

Bayesian age models and stacks: Combining age inferences from radiocarbon and benthic δ18O stratigraphic alignment

Supplement

Model code and software

Viewed

Viewed (geographical distribution)

Cited

2 citations as recorded by crossref.

Bayesian age models and stacks: Combining age inferences from radiocarbon and benthic δ¹⁸O stratigraphic alignment