Bayesian analysis of 210Pb dating
Marco A Aquino-López, Maarten Blaauw, J Andrés Christen, Nicole K. Sanderson
BBayesian analysis of
P b dating
Marco A Aquino-L´opez ∗† Maarten Blaauw ∗ J Andr´es Christen ‡ Nicole K. Sanderson § October 10, 2017
Abstract
In many studies of environmental change of the past few centuries,
P b dating isused to obtain chronologies for sedimentary sequences. One of the most commonly usedapproaches to estimate the ages of depths in a sequence is to assume a constant rate ofsupply (CRS) or influx of ‘unsupported’
P b from the atmosphere, together with aconstant or varying amount of ‘supported’
P b . Current
P b dating models do notuse a proper statistical framework and thus provide poor estimates of errors. Here wedevelop a new model for
P b dating, where both ages and values of supported andunsupported
P b form part of the parameters. We apply our model to a case studyfrom Canada as well as to some simulated examples. Our model can extend beyond thecurrent CRS approach, deal with asymmetric errors and mix
P b with other types ofdating, thus obtaining more robust, realistic and statistically better defined estimates.
Keywords:
P b dating, Chronology, Bayesian Analysis, MCMC, Sediment core. ∗ School of Natural and Built Environment, Queen’s University Belfast, Belfast, BT7-1NN. email: [email protected] , [email protected] † Corresponding author. ‡ Centro de Investigaci´on en Matem´aticas (CIMAT), Jalisco s/n, Valenciana, 36023 Guanajuato, GT,Mexico. email: [email protected] § College of Life and Environmental Sciences, University of Exeter, Exeter, EX4-4QJ. email:
[email protected] a r X i v : . [ s t a t . A P ] O c t Introduction
Radiometric dating is a series of techniques used to date material containing radioactiveelements (Gunten et al., 1995) which decay over time. These techniques use the radioactivedecay equation ( N ( t ) = N e − λt , where N ( t ) is the quantity of a radioactive element left inthe sample at time t , N is the initial quantity, and λ is the element’s decay constant) toinfer the age of the material being investigated. P b (lead-210) is a radioactive isotope which forms part of the U (uranium)series. U (solid) is contained within most rocks and decays into Ra (radium, solid),which later decays into Rn (radon, gas). Since Rn is a gas, a proportion escapes tothe atmosphere where it decays into P b (solid) which is later transported to the earth’ssurface by precipitation.
P b deposited this way is labelled “unsupported” or excess
P b ( P U ). On the other hand, Rn which decays in situ becomes what is labelled “supported” P b ( P S ). By distinguishing between supported and unsupported P b one can determinethe age of the sediment through measuring the
P b at a depth d and compare it to the restof the sediment.The sediments within lakes, oceans and bogs contain layers of biotic and abioticfossils, which can be used as indirect time-series of environmental dynamics as the sedimentsaccumulate over time.. Whereas both unsupported and supported P b decay over time,supported
P b is replenished through decay from radon contained within the sediment.That is why the concentration of supported
P b remains at equilibrium while that of theunsupported
P b decreases and eventually reaches zero. This is the basis for
P b dating.Given its relatively short half-life of 22.3 years,
P b has been used to date recent ( < P b data and a dating model created within a statistical framework, with the objective ofproviding more reliable measures of uncertainty.2
Modelling of
P b data
As outlined above, within sediment
P b is naturally formed from two sources – from sur-rounding sediment and rocks containing U (supported), and from the atmosphere through Rn (unsupported). Modelling these two sources of P b is crucial to the developmentof age-depth models. Since supported and unsupported
P b are indistinguishable fromeach other, in order to model both sources, we have to make assumptions depending on themeasurement techniques used. Measurements of
P b can be obtained by alpha or gammaspectrometry. The latter technique also provides estimates of other isotopes such as Ra ,which can be used as a proxy of the supported P b in a sample (Krishnaswamy et al.,1971).
P b
If gamma spectrometry is used, supported
P b can be assumed to be equal to the con-centrations of Ra . When the sediments are analysed using alpha spectrometry, Ra measurements are not available and estimates of the supported activity can only be ob-tained by analysing sediment which reached background (samples which no longer containunsupported P b ). When alpha spectrometry is used, a constant supported
P b is as-sumed. These two different ways of inferring the supported activity can be formalised bythe following equations: P Ti = P Si + P Ui , (1) P Ti = P S + P Ui , (2)where P Ti is the total P b , P Ui is the unsupported P b and P Si is the supported P b in sample i . Depending on the site and availability of measuring techniques, one of theseequations can be used to differentiate supported from unsupported P b . P b
In order to model the unsupported
P b , some assumptions have to be made regardingthe precipitation of this material from the atmosphere. A reasonable assumption for this3henomenon is the constant flux or rate of supply (Appleby and Oldfield, 1978), whichimplies that for fixed periods of time the same amount of
P b is supplied to the site.Following Appleby and Oldfield (1978), the assumption of a constant rate of supplyimplies that the initial concentration of
P b at depth x (which is linked to age by a function t ( x )), P U ( t ( x )), weighed by the dry mass sedimentation rate r ( t ( x )), is constant throughoutthe core: P U ( t ( x )) r ( t ( x )) = Φ (3)where Φ is a constant. The dry sedimentation rate is the speed at which the sedimentaccumulates, weighed by the sediment’s density at such depth, i.e. r ( t ( x )) = ρ ( x ) dx ( t ) dt , (4)where ρ ( x ) is defined as the density of the sediment at depth x and dx ( t ) dt is the rate at whichthe core accumulates with respect to time. Considering that the relationship between depthand time is expressed by the function t ( x ), then x ( t ) is the inverse function of time, andsince t ( x ) is a one-to-one function r ( t ( x )) = ρ ( x ) (cid:20) dt ( x ) dx (cid:21) − . (5)Since P b is a radioactive isotope it follows from the radioactive decay equation that P U ( x ) = P U ( t ( x )) e − λt ( x ) , (6)where P U ( x ) is the concentration of unsupported P b found at depth x , and λ is the P b half-life. Using equations (3), (5), and (6) the following relationship is obtained, ρ ( x ) P U ( x ) = dt ( x ) dx Φ e − λt ( x ) . (7)Considering that P b is measured over a slice or section of the sediment, this relationshiphas to be integrated over such section to be related to the corresponding measurement, thatis, A U ( a,b ) = (cid:82) ba ρ ( z ) P U ( z ) dz (8)= (cid:82) ba Φ e − λt ( z ) dt ( z ) dz dz = (cid:82) t ( b ) t ( a ) Φ e − λy dy, (9)4here ( a, b ) are the lower and upper depths of the sample respectively and A U ( a,b ) is theactivity in section ( a, b ). Equation (9) provides a link between the age-depth function t ( x )and the unsupported activity in a given section. This is the primary tool to construct anage-depth model based on a constant rate of supply. The Constant Rate of Supply (CRS) (Goldberg, 1963; Robbins, 1978; Appleby and Oldfield,1978) model is the most commonly used
P b dating model. It uses the constant rate ofsupply assumption presented in section 2, and the following equations to obtain a chronology: A U ( x ) = (cid:82) ∞ x ρ ( z ) P U ( z ) dz, (10) A U (0) = (cid:82) ∞ ρ ( z ) P U ( z ) dz, (11) t ( x ) = λ ln (cid:16) A (0) A ( x ) (cid:17) , (12)where A U ( x ) is the remaining unsupported activity below x , and A U (0) is the unsupportedactivity in the whole core. The CRS model can be summarized by equation (12) and fromits term A U (0) one can deduce that this model depends strongly on measuring activitythroughout the whole core. The effect of wrongly estimating this variable is described inAppleby (1998). If the activity cannot be measured throughout the entire core, interpolationis suggested (Appleby, 2001). Moreover, if the lowest sample has not reached background,and thus still contains unsupported P b , extrapolation is suggested.Because this model is based only on the unsupported activity, precise estimates ofsupported
P b are necessary in order to obtain reliable estimates of the unsupported
P b .Depending on the equipment used to obtain the
P b concentrations, and on the modelused to distinguish supported from unsupported
P b , this could be problematic. Wronglyestimating this variable will directly impact the estimate of A (0) which will in consequenceaffect the resulting chronology. 5 .1 Example To show the results of the current approach and later compare them to ours, data obtainedfrom a site in Havre-St-Pierre, Quebec, Canada will be used. The core (HP1C) was obtainedin July 2012 and was analysed using alpha spectrometry at Exeter University, UK . Table 1contains the data from core HP1C. As previously mentioned, alpha spectrometry does notprovide estimates of Ra as is the case for beta spectrometry, but instead, contrary to thelatter, it can measure far smaller quantities of P b . To date this core, the CRS model wascalculated using the recommendations in Appleby (2001).One of the first steps to apply the CRS model is to identify the supported
P b . Forthis purpose the last 4 samples were averaged to obtain an estimate of 8 . Bqkg and a standarddeviation of 1 .
01 for the supported activity. This value was subtracted from the total
P b for each sample, to obtain estimates of unsupported activity. Following Appleby (2001) onecan obtain the dating shown in Figure 1. This methodology requires very strong assumptionsregarding independence, given the fact that it uses accumulated activity as the primary toolfor inference. We now introduce our formal statistical approach for
P b dating, to solvethis and several other issues inherent in the usual CRS technique just described.
P b dating
Let the concentration of
P b in a sample taken from section ( x i − δ, x i ) be a random variable p i . In this case, it is important to clarify that this is not the cumulative concentration oractivity, from the surface to depth x i , but rather it is the concentration found from depth x i − δ to x i where δ is the sample’s thickness. Each concentration of P b ( p i ) is measuredseparately and therefore it is safe to assume that each sample is independent of the othermeasurements and is normally distributed with mean the true concentration P i , and varianceas reported by the laboratory: p i | P i ∼ N (cid:0) P i , σ i (cid:1) . (13)As mentioned above, the supported P b is assumed to be in equilibrium throughout6
Depth (cm) A ge ( y ea r s ) CRS
HP1C
Figure 1: Dating of HP1C obtained by the CRS model (Appleby, 2001) showing the meanand 95% confidence intervals.the core, which means that it remains constant through all depths. If Ra measurementsare available, a supported P b value per sample can easily be included by letting P Si bedifferent at each depth. It is important to note that this will greatly increase the number ofparameters, and should only be used when the hypothesis of a constant supported concen-tration has been shown to be unreasonable. If a constant supported P b is valid, then wecan use equation (2) to infer the supported
P b .Now, assuming a constant rate of supply, as described in section 2, for the unsupported7 P b , the activity in sample i can be obtained as follows A Ui = (cid:90) x i x i − δ ρ i ( z ) P Ui ( z ) dz = (cid:90) t ( x i ) t ( x i − δ ) Φ e − λτ dτ = Φ λ (cid:0) e − λt ( x i − δ ) − e − λt ( x i ) (cid:1) (14)Since the supported P b can be assumed to be constant, the supported activity ofsample i is A Si = (cid:90) x i x i − δ ρ ( z ) P S ( z ) dz = P S ρ i . (15)By defining y i = P Ti ρ i y i | P S , Φ , ¯ t ∼ N (cid:18) A Si + Φ λ (cid:0) e − λt ( x i − δ ) − e − λt ( x i ) (cid:1) , ( σ i ρ i ) (cid:19) . (16)It is important to note that the activity at each sample contains not only the information ofages but also of the supported P b ( P S ) and the initial supply of unsupported P b (Φ)throughout the core.To implement a Bayesian approach, prior distributions for each parameter have to bedefined. Appleby (2001) suggested that the supply of unsupported
P b has a global meanof 50
Bqm yr . This can be used as prior information to obtain a prior distribution for Φ. BecauseΦ is always positive, a gamma distribution can be considered and with this information wecan define Φ ∼ Gamma ( a Φ , b Φ ). On the other hand, since supported P b ( P S ) varies muchfrom site to site, defining an informative prior distribution for P S is in general not possible.As a consequence, data for the supported P b ( P S ) is necessary ( y S , y S , ..., y Sn s ). These datacan come from two different sources; Ra estimates or P b measurements which reachedbackground. A prior distribution for the P S (supported P b ) associated with these datais necessary. Little is known about this parameter prior to obtaining the data. We haveseen cores ranging from nearly 0 up to almost 50
Bq/kg of supported
P b . With thisinformation, a gamma distribution with a mean of 20
Bq/kg and shape parameter of α S = 28ould allow the data to contribute more to the posterior value of P S . Lastly, to define aprior distribution for the ages an age-depth function has to be defined. Since sediment cores can extend back thousands of years,
P b is not the only techniqueused to date them. C (radiocarbon) is a common way to obtain age estimates for organicmaterial. The radiocarbon community has built sophisticated chronology models, whichrely on equally sophisticated age-depth functions, with the objective of reducing and betterrepresenting the uncertainty of the resulting chronology. Because we want our approach tohave the flexibility to incorporate other dating information such as radiocarbon, we decidedto incorporate a well-established age-depth function. Bacon (Blaauw and Christen, 2011), which is one of the most popular chronologymodels for C dating, assumes linear accumulation rates over segments of equal length. Byusing this same structure, age-depth models based on multiple isotopes could potentiallybe obtained. This is the age-depth model we are going to use and we discuss the generalconstruction of the Bacon age-depth function next (see Blaauw and Christen, 2011, fordetails). The age-depth function is considered as linear over sections of equal length, causingdepths to be divided into sections of equal length c < c < c < ... < c K , noting that c = 0. Between these sections linear accumulation is assumed, so for section c i < d < c i +1 the model can be expressed as G ( d, m ) = i (cid:88) j =1 m j ∆ c + m i +1 ( d − c i ) , (17)where m = ( m , m , ..., m k ) are the slopes of each linear extrapolation, and ∆ c = c i − c i +1 isthe length of each section.With this structure a gamma autoregressive model is proposed for the accumulationrate of each section, m j = ωm j +1 + (1 − ω ) α j where α j ∼ Gamma ( a α , b α ) and ω ∈ [0 ,
1] is amemory parameter which is distributed prior to the data as ω ∼ Beta ( a ω , b ω ).Using the above age-depth function, and (16), the log-likelihood for the model takes9he form (cid:96) (¯ y, ¯ y S | m, ω, Φ , P S ) ∝ − (cid:80) ni =1 ( y i − ( A Si + Φ λ ( e − λG ( xi − ,m ) − e − λG ( xi,m ) ))) σ i − (cid:80) n s j =1 ( y Sj − P S )2 σ j . (18)This likelihood contains all the parameters needed by our approach. Using the priordistributions previously mentioned, a posterior distribution f ( m, ω, Φ , P S | ¯ y, ¯ y S ) is defined,from which we intend to Monte Carlo sample the model parameters using MCMC. To allowfor faster convergence of the MCMC, a limit to the chronology is considered. This chronologylimit is inspired by the P b dating horizon, which is the age at which
P b samples lackany measurable unsupported
P b . The
P b dating horizon was described by Appleby (1998) to be 100 - 150 years, basedon the available knowledge and measurement techniques at the time. The dating horizonof a given core is affected by different factors. The first of them is the equipment used tomeasure the samples. If certain equipment has higher precision than another, it will be ableto distinguish unsupported from supported
P b down to deeper samples and thus provideages further back in time. The other factor that affects the dating horizon is the quantity ofinitial unsupported
P b , which is directly affected by the rate of supply (Φ). When there isa larger initial unsupported
P b it will take longer for the unsupported
P b in a sampleto become indistinguishable from the supported
P b .We therefore decided to set a dynamic chronology limit for our method. This limit ( t l )will be determined by two factors – the rate of supply of P b to the site (Φ) and the errorrelated to the equipment used to measure the samples. For example, lets assume that theequipment used to calculate the concentration of
P b in a sample has a minimum error of0 . Bqkg . Now, assuming that the sample comes from a bog with a density ranging between . . gcm (Chambers et al., 2011), then once the unsupported activity in a sample reaches A l (cid:39) . Bqm , it would become indistinguishable from the supported activity. This informationcould help us to calculate the dynamic age limit. By using equation (14) we have A l = (cid:82) t l +1 t l Φ e − λτ dτ = Φ e − λt l − e − λ λ , (19)10here A l is the minimum distinguishable unsupported activity in a sample related to theequipment’s error, Φ is the supply of P b to the site and λ = 0 . − e − λ λ = 0 . t l = λ log (cid:16) . A l (cid:17) (cid:39) λ log (cid:16) Φ A l (cid:17) . (20)It is important to note that this limit depends on the error of the equipment and onthe origin of the samples, which are factors known prior to obtaining the data. Moreover, Φis a variable of the model. This will allow the model to limit the chronology given Φ. Blaauw and Christen (2011) propose the use of a self-adjusting MCMC algorithm, knownas t-walk (Christen and Fox, 2010), which will facilitate the use of these techniques to non-statisticians. The t-walk algorithm requires two initial points for all parameters (Φ , P S , w, α )and the negative of the log posterior function which is called the energy function, U (Φ , P S , w, α | ¯ y, ¯ y S ) = − log f (cid:0) Φ , P S , w, α | ¯ y, ¯ y S (cid:1) . (21)A program (in python 2.7) called Plum is used to implement this approach and to obtaina sample from the posterior distribution.
Plum has been tested on peat and lake sedimentcores, as well as on simulated data, providing reasonable results with no tuning of theparameters needed; examples of these results can be seen in sections 5 and 6. The consistencyof these results, with minimal user input, show how the t-walk (Christen and Fox, 2010) wasa suitable choice for this implementation.
To implement our approach to the HP1C data presented in section 3.1,
Plum was pro-grammed to use the last 4 samples from Table 1 as estimates of the supported activity, usingthe rest of the samples to establish the chronology. Figure 2 shows the results of the CRSmodel in red and our approach in black and grey. From this comparison we can observe11hat both models agree with each other down to a depth of 25 cm, at which point the CRSmodel continues at a similar slope unlike our approach which provides younger estimates.This uninterrupted growth of the CRS model can be explained by its age function which isa logarithmic approximation, invariably tends to infinity as unsupported
P b reaches 0 .Even with these discrepancies both models have overlapping confidence intervals, with ourapproach providing a more precise chronology in the topmost part and a more conservativeestimate for the deepest part of the core.This example shows the potential of our approach in a ‘well-behaved’ real-world casestudy, but unfortunately we cannot observe the precision of this approach when confrontedwith more challenging data sets, such at those which did not reach background and/or withmissing data. For this purpose, several simulations were created where we know the ‘true’chronology and can observe how our approach behaves in more challenging circumstances.12
Depth (cm) A ge ( y ea r s ) PlumCRS
HP1C
130 135 140 145 . . . . Supply D en s i t y l )( 5 6 7 8 9 10 11 12 . . . . Supported Activity
BqKg D en s i t y l )( Figure 2: Comparison between the CRS (Appleby, 2001) and our model using data fromHP1C. Blue curve and shadow indicate CRS mean and its corresponding 95% range. Dashedblack curves indicate mean and 95% confidence interval for our model. Grey lines are sim-ulations from
Plum . The top curves represent estimates of the supply of unsupported
P b (Φ) and supported
P b ( P S ) using the CRS model (red; dot shows the mean, parenthesesshow the standard deviation) and Plum (black curve).13
Simulated Example
To obtain simulated data, a constant supply of
P b was defined as Φ = 150
Bqkg , and byusing the constant rate of supply assumption from equation (3) we have P ( x ) r ( x ) = 150.At this point, we can define ρ ( x ) to obtain r ( x ) by using equation (5) and the age function t ( x ) = x / x/ ρ ( x ) = 1 . − .
05 cos (cid:16) x π (cid:17) (22) P ( x ) = 150( x + ) ρ ( x ) . (23)Using these functions, simulated samples at any given depth can be obtained byintegrating each function between the top and bottom depths of the sample. Lastly, tosimulate supported P b a constant value was added to the simulations such that P i = P S + (cid:82) ba P U ( x ) dx , where a and b are the top and bottom depths of the sample. For thissimulation we set the supported P b to P S = 20. To replicate the measurement errorsrelated to the concentration of P b , white noise was added such that P i + (cid:15) where P i isthe concentration found in sample i and (cid:15) ∼ N (0 , σ i ). This exercise provided us with thedataset in Table 2. We use this simulated data set to test the precision of our approachin various circumstances. For this purpose, the last three sample points were designated asestimates of the supported P b .The best scenario for
P b age-depth models is when every core section is measured,from the surface to where background is reached. In this scenario any approach should reachthe best results, thus providing the complete information about the decay of unsupported
P b . This scenario can be simulated using the complete data set from Table 2. Figure3 shows the comparison between the chronology obtained by our model and that of theCRS Appleby (2001) alongside the real age function, and how both models include thetrue chronology in their 95% intervals. By applying our approach to this scenario, weobtained a very accurate chronology by taking the mean of the MCMC simulations. Thisshows, unsurprisingly, that our model behaves quite well in the best-case scenario. On theother hand, the CRS model provides a shorter chronology, because some samples had to bediscarded from the chronology. This is a direct result from the logarithmic approximationmentioned in section 5. In this particular case, the two bottommost samples had to be14iscarded since the last sample was smaller than the mean of the three samples used tocalculate the supported activity. On the other hand, CRS estimates younger ages for thisexample, which can be a result of the underestimated supported
P b as can be observed infigure 3. Another feature of the CRS worth mentioning is the rapid growth of the chronologyin the last sample. As previously mentioned, this rapid increase can be attributed to thelogarithmic approximation the CRS uses.
Simulated data
Depth (cm) A ge ( y ea r s ) True age modelPlumCRS
10 15 20 25 30 35 . . . . Supported Activity
BqKg D en s i t y )( Figure 3: Comparison between the CRS (Appleby, 2001) and our model using simulated data.Blue curve and shadow indicate CRS mean and its corresponding 95% range. Dashed blackcurves indicate mean and 95% confidence interval for our model. Grey lines are simulationsfrom
Plum . Red curve is the true age-depth model. The curve in the right represent estimatesof the supported
P b ( P S ) using the CRS model (red; dot shows the mean, parenthesesshow the standard deviation) and Plum (black curve). True supported
P b ( P S ) is markedby a blue line.The following scenarios deal with the behaviour of our model in circumstances wherethere is not complete dating information. Even if we attempt to use the CRS model toprovide age estimate in these scenarios, it does this by interpolating and extrapolating inthe sections where there is missing data. Applying the CRS model to these simulations wouldrequire us to take several additional heuristic decisions with large potential impacts on the15hronology (e.g., exponential or linear extrapolation to beyond and/or between the datedlevels etc., see Sanchez-Cabeza and Ruiz-Fern´andez (2012)). Such comparisons lie outsidethe scope of the present work but will be explored in a future study and consequently forthe next examples we only study the performance of the Plum chronology.Sometimes researchers do not have the funds to obtain a full, continuously measureddataset for the chronology that they want to build. When this is the case, only certainstrategically placed samples are measured. To simulate this scenario, only the data at odddepths was used to obtain the chronology. Figure 4 shows the results from this experiment.The accuracy of the model did not change as it still gives an accurate estimate of the trueage model, and the precision was not greatly affected even though only half of the availabledata was used to calculate this chronology.A common problem in P b dating is not reaching background. To observe thebehaviour of our model in these circumstances, the bottommost seven data points wereremoved leaving us with a dataset which has not reached background. Figure 4 presentsthe resulting chronology compared to the true age function. The chronology seems to beaccurate down to a depth of 16 cm, from which point it seems to provide older estimates.On the other hand, the model encloses the true chronology at all times even for the olderages. The last scenario to which our approach was tested is missing topmost sediment. Forthis example, the data points with a depth of two to ten cm depth were removed leaving uswith a data set with topmost missing data. Figure 4 shows the results of this experiment.Even with a third of consecutive missing data, the model is able to accurately reconstructthe true age function.Our approach behaves well in every tested scenario, as its accuracy is not greatlyaffected by any of the different scenarios we introduced.
The approach developed here presents a more robust methodology to deal with
P b data.The advantage comes from a more realistic measure of uncertainties, since the ages are16
Missing topmost sediment
Depth (cm) A ge ( y ea r s ) True age modelPlum
Strategic missing data
Depth (cm) A ge ( y ea r s ) True age modelPlum
Background not reached
Depth (cm) A ge ( y ea r s ) True age modelPlum
Figure 4: Bayesian analysis of simulated
P b data using odd depths in the top-left, usingsamples with depths 1-20 in the top-right and using the samples with depths 1 and 11-27in the bottom-centre. The red line represents the true age-depth function, grey lines aresimulations from
Plum ; dashed lines represent the 95% interval and mean.parameters which are inferred in the process. Moreover, dealing with missing data, whichis a common problem when dealing with
P b dating, becomes easier because our modeldoes not need the whole core to be measured to obtain accurate results. Also, since the CRSmodel relies on a ratio, that approach requires removal of the bottommost measurement.Since our methodology does not rely on a ratio, all the samples provide information to the17hronology, making longer chronologies possible.Because of the integration of the supported
P b into the model a posterior distribu-tion of this variable can be obtained, as well as for ages at any given depth (not just thosewith
P b measurements) and the supply of
P b to the site. Figure 2 shows the posteriordistributions of the supported
P b and the supply of unsupported
P b . These posteriordistributions provide more realistic estimates of the uncertainty of these variables, whichmay be used for other studies where the main focus is not the chronology but other aspectsof the site.Another advantage of this methodology is the fact that since the model operates withina Bayesian framework, incorporating extra information is possible without having to ‘double-model’ by using previously modelled ages within an age-depth model. This informationcould come in the form of other radiometric ages, such as radiocarbon determinations. Sincemeasurements of radiocarbon and
P b , given the age, are independent, the overall likelihoodwould consist of two parts; the likelihood from
P b and from C . Therefore, L (Θ) = L P b (Θ) L C (Θ) . (24)Considering that the only link between both data is t ( x ), by using the same age-depth function such as that from equation (17), a chronology with both sources of data ispossible. This becomes very important because the calibration curve (Reimer et al., 2013),which is used to correct the radiocarbon ages, is non-linear for the most recent few centuries,causing problems with interpreting radiocarbon ages. This period is partly covered by P b .By combining these two methodologies, more robust chronologies can be obtained for thisimportant period in human and environmental history.
References
Appleby, P. G. (1998), “Dating recent sediments by Pb-210: Problems and solutions,”
Proc.2nd NKS/EKO-1 Seminar, Helsinki, 2-4 April 1997, STUK, Helsinki , pp. 7–24.Appleby, P. G. (2001), “Chronostratigraphic Techniques in Recent Sediments,”
Tracking nvironmental Change Using Lake Sediments: Basin Analysis, Coring, and ChronologicalTechniques , pp. ”171–203”.Appleby, P., and Oldfield, F. (1978), “The calculation of lead-210 dates assuming a constantrate of supply of unsupported 210Pb to the sediment,” Catena , 5(1), 1–8.Blaauw, M., and Christen, J. A. (2011), “Flexible paleoclimate age-depth models using anautoregressive gamma process,”
Bayesian Analysis , 6(3), 457–474.Chambers, F. M., Beilman, D. W., and Yu, Z. (2011), “Methods for determining peathumification and for quantifying peat bulk density , organic matter and carbon contentfor palaeostudies of climate and peatland carbon dynamics,”
Mires and Peat , 7, 1–10.
URL: http://mires-and-peat.net/pages/volumes/map07/map0707.php
Christen, J. A., and Fox, C. (2010), “A general purpose sampling algorithm for continuousdistributions (the t-walk),”
Bayesian Anal. , 5(2), 263–281.
URL: http://dx.doi.org/10.1214/10-BA60
Goldberg, E. D. (1963), “Geochronology with Pb-210,”
Radioactive Dating , pp. 121–131.Gunten, B. H. R. V., Institut, P. S., and Psi, C.-V. (1995), “Radioactivity: A Tool to Explorethe Past I . Introduction,”
Radiochimica Acta , 70/71, 305–316.Krishnaswamy, S., Lal, D., Martin, J., and Meybeck, M. (1971), “Geochronology of lakesediments,”
Earth and Planetary Science Letters , 11(1-5), 407–414.Reimer, P. J., Bard, E., Bayliss, A., Beck, J. W., Blackwell, P. G., Bronk, M., Grootes, P. M.,Guilderson, T. P., Haflidason, H., Hajdas, Hatte, C., Heaton, T. J., Hoffman, D. L., Hogg,A. G., Hughhen, K. A., Kaiser, J. F., Kromer, B., Manning, S. W., Niu, M., Reimer,R. W., Richards, D. A., Scott, E. M., Southon, J. R., Staff, R. A., Turney, C. S. M., andvan der Plicht, J. (2013), “IntCal 13 and Marine 13 radiocarbon age calibration curves050,000 years cal BP,”
Radiocarbon , 55, 1869–1887.Robbins, J. (1978), “Geochemical and geophysical applications of radioactive lead,”
Thebiogeochemistry of lead in the environment , pp. 285–393.19anchez-Cabeza, J. A., and Ruiz-Fern´andez, A. C. (2012), “210Pb sediment radiochronology:An integrated formulation and classification of dating models,”
Geochimica et Cosmochim-ica Acta , 82, 183–200.Sanchez-Cabeza, J. A., Ruiz-Fern´andez, A. C., Ontiveros-Cuadras, J. F., P´erez Bernal,L. H., and Olid, C. (2014), “Monte Carlo uncertainty calculation of 210Pb chronologiesand accumulation rates of sediments and peat bogs,”
Quaternary Geochronology , 23, 80–93. 20
Data
Depth
P b σ ( P b ) Density ( ρ ) Depth P b σ ( P b ) Density ( ρ ) cm Bq/kg σ ( Bq/kg ) g/cm cm Bq/kg σ ( Bq/kg ) g/cm .
730 11 .
900 0 .
045 18 279 .
320 11 .
140 0 . .
390 15 .
080 0 .
047 19 243 .
820 9 .
940 0 . .
240 17 .
110 0 .
051 20 246 .
750 9 .
170 0 . .
640 14 .
430 0 .
049 21 351 .
680 13 .
100 0 . .
040 16 .
440 0 .
049 22 281 .
280 11 .
380 0 . .
970 16 .
750 0 .
051 23 235 .
300 12 .
720 0 . .
120 16 .
780 0 .
050 24 192 .
820 7 .
240 0 . .
880 15 .
200 0 .
047 25 94 .
280 4 .
740 0 . .
580 14 .
830 0 .
048 26 50 .
550 3 .
410 0 . .
750 16 .
210 0 .
049 27 36 .
080 2 .
260 0 . .
310 13 .
030 0 .
052 28 28 .
710 2 .
100 0 . .
770 15 .
220 0 .
047 29 24 .
680 1 .
760 0 . .
740 13 .
450 0 .
051 35 11 .
040 1 .
270 0 . .
410 10 .
020 0 .
050 40 6 .
240 1 .
010 0 . .
580 11 .
620 0 .
053 45 10 .
150 1 .
310 0 . .
170 9 .
760 0 .
048 50 7 .
960 1 .
600 0 . .
740 12 .
950 0 . P b ( P T ) σ ( P b ) Density ( ρ ) Depth P b ( P T ) σ ( P b ) Density ( ρ ) cm Bq/kg σ ( Bq/kg ) g/cm cm Bq/kg σ ( Bq/kg ) g/cm .
897 10 0 .
145 16 80 .
845 7 0 . .
761 9 0 .
145 17 64 .
024 7 0 . .
507 9 0 .
145 18 48 .
792 7 0 . .
669 9 0 .
145 19 54 .
076 7 0 . .
026 9 0 .
146 20 37 .
109 7 0 . .
949 9 0 .
146 21 36 .
640 7 0 . .
226 9 0 .
146 22 28 .
602 7 0 . .
736 9 0 .
146 23 22 .
180 6 0 . .
598 8 0 .
147 24 29 .
342 6 0 . .
080 8 0 .
147 25 28 .
723 6 0 . .
818 8 0 .
148 26 26 .
123 6 0 . .
937 8 0 .
148 27 17 .
803 6 0 . .
476 8 0 .
149 28 23 .
349 6 0 . .
268 8 0 .
149 29 13 .
607 6 0 . .
069 8 0 .
150 30 16 .