- Email: [email protected]
A new probabilistic technique to build an age model for complex stratigraphic sequences
Quaternary Geochronology 22 (2014) 65e71
Contents lists available at ScienceDirect
Quaternary Geochronology journal homepage: www.elsevier.com/locate/quageo
Research paper
A new probabilistic technique to build an age model for complex stratigraphic sequences Martin H. Trauth* Institut für Erd- und Umweltwissenschaften, Universität Potsdam, 14476 Potsdam, Germany
a r t i c l e i n f o
a b s t r a c t
Article history: Received 25 June 2013 Received in revised form 18 February 2014 Accepted 4 March 2014 Available online 20 March 2014
The age models of fluvio-lacustrine sedimentary sequences are often subject of discussions in paleoclimate research. The techniques employed to build an age model are very diverse, ranging from visual or intuitive estimation of the age-depth relationship over linear or spline interpolations between age control points to sophisticated Bayesian techniques also taking into account the most likely deposition times of the type of sediment within the sequence. All these methods, however, fail in detecting abrupt variations in sedimentation rates, including the possibility of episodes of no deposition (hiatus), which is the strength of the method presented in this work. The new technique simply compares the deposition time of equally thick sediment slices from the differences of subsequent radiometric age dates and the unit deposition times of the various sediment types. The percentage overlap of the distributions of these two sources of information, together with the evidence from the sedimentary record, helps to build an age model of complex sequences including abrupt variations in the rate of deposition including one or many hiatuses. Ó 2014 Elsevier B.V. All rights reserved.
Keywords: Age-depth modelling Stratigraphy Sedimentation rate Hiatus MATLAB
1. Introduction A precise age-depth model is an important prerequisite for the interpretation of high-resolution paleoclimate records. Unfortunately, these methods, although often published with downloadable computer code, guidelines and work-through examples, are hardly used in paleoclimate research (Blaauw, 2010). Many paleoclimate studies do not provide any information on the algorithm used in calculating the age model, if they ever use one, or they do not provide any information about the ambiguity and uncertainties of the result. The available age-depth modelling techniques used range from (1) visual or intuitive estimation of the age-depth relationship (e.g. Brown and Feibel, 1991; Trauth et al., 2001; Behrensmeyer et al., 2002), (2) linear, polynomial or spline interpolation or regression between the radiometric age dates (e.g. Maher, 1972; Blaauw and Heegaard, 2012), (3) calibrating the stratigraphy to insolation or orbital target curves (e.g. Partridge et al., 1997; Joordens et al., 2011), (4) forward modelling of facies variations in cyclic sections (e.g. Kominz and Bond, 1990, 1992), (5) Monte-Carlo modelling of age distributions along sediment cores (Hercman and Pawlak, 2012) to (6) Bayesian age-depth models (e.g.
* Tel.: þ49 331 977 5810; fax: þ49 331 977 5700. E-mail addresses: [email protected], [email protected]. http://dx.doi.org/10.1016/j.quageo.2014.03.001 1871-1014/Ó 2014 Elsevier B.V. All rights reserved.
Steier and Rom, 2000; Steier et al., 2001; Blaauw and Christen, 2011; Blaauw and Heegaard, 2012). The weakness of many of these techniques is mainly due to a lack of a suitable algorithm for consideration of abrupt changes in the sedimentation rate within stratigraphic sections or cores, including even the possibility of the complete lack of sedimentary layers at certain times (hiatus). Building an age model in these situations is complicated by the fact that the sedimentation rates depend strongly on the observed time resolution, with a very strong negative relationship between expected sedimentation rate and averaging time (e.g. Sadler, 1981, 1999) (Fig. 1). Such changes, however, are very common in stratigraphic sequences, for instance in the tectonically-active sedimentary basins of Eastern Africa (e.g. Brown and Feibel, 1991; Trauth et al., 2001; Behrensmeyer et al., 2002). In these basins, the sedimentation rates range from less than 0.1 m kyr1 for diatomite (a sediment composed of the skeletons of silica algae, see Table 1), and 0.1e10 m kyr1 for clastic sediments in lakes (Einsele, 2000; Hinderer and Einsele, 2001), to more than >10 m/kyr for sands (Einsele, 2000), and several meters of volcanic air fall deposits per within a couple of hours, followed by a longer time of no deposition (Fisher and Schmincke, 1984). Whereas the sedimentation rates of most sediment types are well known for lakes in Eastern Africa, the significance of hiatuses are subject passionate discussions (Trauth et al., 2005; Trauth and Maslin, 2009; Owen et al., 2009).
66
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
Fig. 1. Possible age models for three layers of clay/silt and diatomite, constrained by two radiometric age dates, 100 and 400 kyr, and their errors. The difference age models include different sedimentation rates, even with similar types of sediment, as well as one or more hiatuses of different duration. The sedimentation rates of the individual types of sediment are averaged over the entire interval. As an example, the lower clay/silt bed in Age Model 2 has a sedimentation rate of 10 m/150 kyr z 0.07 m kyr1 ¼ 0.07 m kyr1, that of the upper clay/silt bed is 10 m/50 kyr ¼ 0.2 m kyr1 ¼ 0.2 m kyr1, which averages as 20 m/200 kyr ¼ 0.1 m kyr1 ¼ 0.1 m kyr1. The choice of the age model has important implications for the sedimentation rate of the diatomite and hence for the duration of lake episode documented by the diatomite.
The possibility of an abrupt change in the sedimentation rate was even regarded as a weakness of a linear over a spline age model, even though these changes exist but are often excluded or ignored in stratigraphic sections, similar to the possibility of hiatuses. Instead, an age model has to be smooth, often by introducing an arbitrary chosen memory in the sedimentation rates along the core, avoiding extreme variations or extremely low or high sedimentation rates, without abrupt changes in the measures of central tendency and dispersion, and without major gaps (e.g. Bronk Ramsey, 2008, 2009; Blaauw and Christen, 2011). Furthermore, it is assumed that a single sediment type, defined by the Gammadistributed sedimentation rate, occurs through the section, which could be a good model in some lacustrine settings, but is an undue assumption in many sedimentary environments (e.g. Bronk Ramsey, 2008, 2009; Blaauw and Christen, 2011). Here, I present a new technique to determine the age-depth relationship for stratigraphic sections and sediment cores with
extremely fluctuating sedimentation rates, including hiatuses. The technique uses two independent sources of information to build an age model, (1) the time of deposition based on radiometric ages with their Gaussian errors and (2) the typical unit deposition time of the various sediment types (Fig. 2). Although the algorithm has been designed, tested and applied to age dates with Gaussian errors, it can be adapted for other applications quite easily. The agreement, or disagreement, between the two estimates of the time of deposition from two statistically independent sources of information helps to build an age-depth relationship in complex stratigraphic sequences. While this new technique in the case of high-resolution radiocarbon chronologies of relatively young (<104 yrs) and homogeneous sediments is not as straight forward as Bayesian techniques, it is superior in older 40Ar/39Ar dated sedimentary sequences with extremely fluctuating sedimentation rates and high probabilities for the occurrence of hiatuses. 2. The system to be modelled
Table 1 Typical sedimentation rates of diatomite as documented in Quaternary lake sediment sequences from East Africa. Location Gadeb, Ethiopia Ol Njorowa Gorge Ol Njorowa Gorge Gicheru Olorgesailie Gicheru Lake Malawi Kariandusi Lake Manyara Barsemoi
Sedimentation rate 1
0.1 m kyr 0.03e0.25 m kyr1 0.13e0.4 m kyr1 0.5 m kyr1 1.0 m kyr1 1.2 m kyr1 >1.3 m kyr1 1.5 m kyr1 1.9 m kyr1 2.2 (0.3e12) m kyr1
Reference Gasse 1980 Bergner and Trauth, 2004 Trauth et al., 2001 Trauth et al., 2007 Deino and Potts, 1990 Trauth et al., 2007 Owen and Crossley, 1992 Trauth et al., 2007 Holdship, 1976 Deino et al., 2006
Sedimentary sections consist of a layered sequence of materials with different grain size and composition. The rate of deposition of these materials is described by the sedimentation rate (sediment thickness per unit time) or, alternatively, the accumulation rate (solid sediment mass per unit area and time, in kg m2 kyr1) (Einsele, 2000). The sedimentation rate (in m kyr1) is determined from the sedimentary unit of the thickness z, divided by the time difference t between two radiometric age dates below and above the unit. While the terms sedimentation and accumulation rate are usually used correctly in the literature, there is no consensus on the name of the inverse of the sedimentation rate (time per sediment thickness). I prefer the term unit deposition time (in kyr m1), which
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
Fig. 2. (A) Gaussian frequency distribution of the time difference per slice between two age dates with Gaussian errors is determine by radiometric dating. (B) The Gamma distributions of the unit deposition time defined for each of the type of deposit (e.g. tephra, clastic sediments, and diatomite) is determined by the sediment record.
is very close to deposition time used by Goring et al. (2012) but different from accumulation rate used by Blaauw and Christen (2011). The sedimentation rate varies dramatically between different depositional environments, ranging from >104 m kyr1 in continental slopes, alluvial fans and delta plains to <101 m kyr1 in lakes and the deep sea below the carbonate compensation depth (Einsele, 2000). Of course, these rates also vary within the same environment, even at the same place, driven by volcano-tectonic, climate and other factors, especially in tectonically active regions (Einsele, 2000). Unfortunately, the stratigraphic record is fragmentary, i.e. it contains stratigraphic gaps without deposition, often without any information on the exact duration of these episodes of no deposition or hiatuses (Sadler, 1981, 1999; Sommerfield, 2006; Schumer and Jerolmack, 2009). The possibility of fluctuating sedimentation rates and several hiatuses between two radiometric age dates complicates the assessment of the actual time of deposition of a sedimentary sequence. As an example, let us assume a simple sequence of three sediment layers, a unit of clastic sediments overlain by a diatomite bed, which is again overlain by clastic sediments, each of these sediment layers being 10 m thick (Fig. 1). The sedimentary sequence is bounded by two age dates, 100 and 400 kyr old. Fig. 1 shows seven possible age models for the sedimentary sequence, all of them consistent with the plausible sedimentation rates of fine-grained sediments (0.0697e0.1 m kyr1) and diatomite (0.1e10 m kyr1). Age models 3, 4 and 6, however, also have hiatuses with a duration of 40e140 kyr, which complicates the interpretation of the stratigraphic sequence. The hiatuses can occur both at the layer boundaries (age model 6) as well as within layers (age models 3 and 4). How can we find the best (or best-fitting) age model for such sedimentary sequences? In statistics, or in sciences in general, the golden rule is always to use the simplest model, which is explains the data. In mathematical terms, the model with the smallest
67
number of adjustable parameters is considered the best (or bestfitting) as long as the mean least-squared error is below a critical value. As an example, if a simple model with a linear sedimentation rate (two parameters, slope and y-intercept) describes the data well in a least-squared sense, we prefer it over a model with three parameters (cubic polynomial, three parameters). In practice, a piecewise linear or cubic spline interpolation between age control points is used to find the best age model (e.g. Blaauw and Heegaard, 2012). In this context, it has often been argued that cubic spline interpolation between age control points is better than simple linear interpolation because it smoothes out minor errors in age data points (Maher, 1972; Blaauw and Heegaard, 2012). This implies that abrupt changes in sedimentation rates are due to errors in the age dates. On the other hand, in more complicated sequences with varying sediment types, such abrupt changes in the sedimentation rates may actually exist, as shown in Fig. 1, also showing zero sedimentation rate (hiatuses); even negative sedimentation rates (erosion) are possible. In fact, the advantage of the technique presented here is the ability to detect such episodes without deposition in the record without any prior assumptions. Furthermore it allows to include multiple sediment types with different sedimentation rates. Most importantly, there are no a priori assumptions made regarding a memory or coherence in the sedimentation rates in the age model to suppress abruptly changing rates; abruptly changing rates are explicitly expected at discontinuities between layers of different sediment types. In such an age model, a linear interpolation between the age control points could provide the better age model over cubic splines. 3. Description of the algorithm I use a new technique to build an age model for stratigraphic sections and sediment cores (Fig. 1). The basic idea of the new method is to estimate the most likely time of deposition between two radiometric age dates, such as 40Ar/39Ar dates, based on two sources of information (Fig. 2). First, the time difference is determined from the two age dates and their independent Gaussian errors. Second, the time difference is again calculated from the time accumulation rate of the sediment types and their Gammadistributed dispersion. The time difference per slice between two age dates is determined by radiometric dating of volcanic ash layers in the case of 40Ar/39Ar dates of feldspar crystals or any other dating technique (Fig. 2a). The age dates are statistically independent from the type of sediment (or any other deposit) between the dated levels in the stratigraphic section. The errors of the ages are assumed to be independent from each other. Therefore, the time difference between two age dates can be calculated, for instance by using the Gaussian law of error propagation. In practice, calculating the difference of two ages in a stratigraphic section is slightly more difficult than described here. Radiometric age dates such as 40 Ar/39Ar ages of subsequent sediment layers, or more precisely, the one of the intercalated ash layers, are often not significantly different within the error bars (zero age difference) or the estimated age of the stratigraphically higher layer is older than the lower one (age reversal). For simplicity in our experiment, we assume monotonically increasing ages with depth and therefore exclude those ages that are unexpectedly younger or older than expected from the stratigraphic sequence. Of course we have to explain these outliers or design a more complex age model for the section, i.e. considering reworking of certain units in the section. The unit deposition time is defined for each of the type of deposit in the stratigraphic section (Fig. 2b). First, the type of sediment is defined in the outcrop or in the sediment core. We can work with
68
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
any number of sediment types, but at the beginning we should work with a small number. Second, a Gamma model is used for the most likely unit deposition time for each sediment type (e.g. three types of materials: tephra, clastic sediments, diatomite). The Gamma distribution has two parameters, shape and scale. In the experiments described here, a Gamma model gamma(shape, scale) with shape ¼ 2 best describes the frequency distribution of unit deposition times, which is in agreement with the work of Blaauw and Christen (2011). At the lower end of unit deposition times is diatomite with typical sedimentation rates between 0.03 and 12 m kyr1, with a mode at around 1 m kyr1 (Table 1). In our experiments, a Gamma model with gamma(2,1.5) best describes the expected range of unit deposition times. For clastic sediments we found gamma(2,2) as being the best-fitting model, corresponding to the most frequent sedimentation rate of 0.5 m kyr1. The best Gamma model for tephra is gamma(2,8) with a very wide range of possible unit deposition times, taking into account the large variation in possible sedimentation rate between several millimetres of fine-grained ash and thick layers of air fall tephra deposited within hours. The agreement, or disagreement, between the two estimates of the time of deposition from two statistically independent sources of information helps to build an age-depth relationship in complex stratigraphic sequences. If the two estimates of the time of deposition match, the probability of a hiatuses is very low. The time difference between the two estimates is regarded as the duration of the hiatuses. Because both values are linked to their frequency distributions, we can also calculate the error probability for the conclusion “hiatus” and its duration. While this new technique in the case of high-resolution radiocarbon chronologies of relatively young and homogeneous sediments is not as straight forward as Bayesian techniques, it is superior in sedimentary sequences with fluctuating sedimentation rates and high probabilities for the occurrence of hiatuses. Alternatively, if there is no evidence for an hiatus between the age dates, we can multiply the two independent frequency distributions, Gaussian for the age difference and Gamma for the sediment types, to calculate a better estimate of the time of deposition and its frequency distribution. After having determined the distribution of the deposition times along the stratigraphic record, the paleoclimatologist wants to know the best-fitting age model at which one can interpolate the physical, chemical and biological environmental indicators. Much has been written about using the mean, weighted mean, median, mode, or any other more complex measure of central tendency to build the age model from the frequency distribution. I found the argument of Blaauw and Heegaard (2012) very well justified, recommending to use the mode (modal value) of the probability distribution over any other measure of central tendency. According to them, this has the advantage of being the age that is most likely, but it can concentrate emphasis on one peak that happens to be fractionally higher than another in a distribution that is bimodal or even multi-modal, which could certainly be the case in radiocarbon-dated sections rather than in Argon-dated sequences. Similar to them, I use percentiles of the probability distribution to describe the spread or uncertainty used as tolerance, for example, 25% below and above the mode in the example below. 4. Examples 4.1. Synthetic Argon-dated lake-sediment sequence Here, I demonstrate the use of the new technique to build an age model for a synthetic 50 m thick lake sediment sequence (Fig. 3). The Supplementary material of the paper contains a MATLAB script called script.m and the input files input_section.txt and
Fig. 3. Best-fitting age model of a synthetic section with three sediment types (tephra, clastic sediments, and diatomite) and radiometric age dates (in black, with one-sigma error bars). The age model (blue line) was calculated from the modes of the multiplied distributions of the age differences per slice and the unit deposition times of the deposits, displayed together with the 25% confidence bounds below and above the modes (red lines). Based on the initial offset of the age model and the radiometric age dates above 37.4 m height and the relatively low percentage overlap (w7.4%) of the distributions, a 2 kyr long hiatus was introduced at 35 m. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
input_ages.txt as an example to the application of the technique. The MATLAB script contains all necessary steps to build the age model based on the lithological classification and ages of the sedimentary sequence in the input files. The input file input_ section.txt has the format
0 6 12 ð.Þ 38 42 47
6 12 14
1 3 1
42 47 50
1 2 1
where the first and second column include the base and top of each sedimentary unit in meters. The third column includes an identifier for the type of deposit, in which 1 stands for tephra, 2 for clastic sediments and 3 for diatomite. The input file input_ages.txt contains the Argon ages in the format
2:8 398:8 2:0 7:6 394:0 1:5 11 387:4 2:0 ð.Þ 37:4 356:7 3:0 43:2 352:8 1:6 47:8 348:8 1:2 where the first column includes the height of the sample in the section in meters. The second and third column contains the ages before present and their one-sigma errors. After loading the data into the workspace of MATLAB, we define a starting value for the maximum age maxage of the section. The initial choice of the value of maxage is not critical to the result of the experiment and can be modified in a later run of the model to get better results. The section is then divided into equally thick slices of sediment, for which the deposition time is determined from the radiometric ages and the unit deposition time of the sediment types. The choice of the thickness slices defines the observed time resolution and
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
therefore the amplitude of possible variations in the unit deposition times. Furthermore, smaller values for the thickness slices produce more accurate results at the expense of higher computation time. In the example, I used slices ¼ 0.1 for 0.1 m thick slices to get a good result. The following set of commands in the MATLAB script includes the interpolation of the sedimentary data from input_section.txt upon an evenly spaced grid of the 0.1 m slices. The first estimate of the deposition time is based on the unit deposition time of the sediment types in the section. We first design the Gamma models for the three sediment types in the section as classified in input_section.txt. Similar to the Bayesian techniques by Blaauw and Christen (2011), we can start with an initial set of shape and scale factors and then we iteratively improve the factors to get better results. In the synthetic experiment, I used a shape factor of 2 for all sediment types and a scale factor of 0.5, 0.8, and 2.0 for the three sediment types tephra, clastic sediments, and diatomite, respectively, as defined in the array shape ¼ [2 0.5;2 0.8;2 2]. The Gamma models for the sediment types helps to calculate a series of unit deposition times along the section, from which we calculate the modes and the quartiles, as well as the cumulative deposition time of the section. The second estimate of the deposition time comes from the difference of the Argon ages in the file input_ages.txt. After calculating the difference of subsequent ages, we interpolate a series of deposition ages along the section for each sediment slice, calculate the modes and the quartiles, as well as the cumulative deposition time of the section. Since we now have two independent estimates for the deposition time including the normalised Gamma/Gaussian distributions of these estimates, we can calculate the percentage overlap of the normalised distributions. The value of 1eoverlap helps to define intervals, which probably contain one or several hiatuses. If there is no reason to prefer one of the two estimates for the deposition time, then we can multiply the two distributions to obtain an improved estimate of the true deposition time for each sediment slice. In our example, we notice very low values (w7.4%) of the percentage overlap of the two distributions between 33.9 and 36 m. Here, the deposition time from the sediment types is 0.33, whereas the deposition time from the radiometric ages is 0.37 years. In contrast, the deposition time at around 10 m height is in good agreement, both being w0.19 years per 0.1 m thick slice. To compensate for the discrepancy between the two values for the deposition time, we can now install a hiatus between 33.9 and 36 m in the section, up to the full difference between the ages and the expected unit deposition time of the sediment. The exact location of the hiatus between 33.9 and 36 m is unknown, but may be found by looking into the stratigraphic logs of the outcrop. In this synthetic example, I have decided to introduce a 2 kyr long hiatus in the middle of the interval at 35 m, which gives a good age-depth model for the section. 4.2. Argon-dated lake-sediment sequence of Lake Naivasha Here, I use the new technique to build an age model for the lake sediment sequences, ca. 150e60 kyr old, exposed in the Ol Njorowa Gorge south of Lake Naivasha in the Central Kenya Rift (Trauth et al., 2001, 2003) (Fig. 4). Anorthoclase and sanidine phenocryst concentrates from 16 tephra beds were dated by the laser-fusion 40 Ar/39Ar method. The chronology is free of age reversals within the one-sigma error bars of the 40Ar/39Ar ages. The means, however, show in fact two reversals at the top (61.9 m, 60 2 kyr) and in the middle of the section (40.28 m, 92 5 kyr). Since my technique, similar to BACON by Blaauw and Christen (2011), requires
69
Fig. 4. Best-fitting age model of the Late Pleistocene lacustrine section exposed in the Ol Njorowa Gorge, south of Lake Naivasha, Kenya. The w60 m thick sequence consists of three types of deposits, (1) tephra, (2) clastic sediments, and (3) diatomite, based on the stratigraphic column of the section published by Trauth et al. (2001, 2003). Based on the percentage overlap of the distributions of the age differences per slice and the unit deposition times of the deposits, also considering the presence of the 10 m thick Kedong Valley Tuff (107 kyr old) in the middle of the section, two hiatuses, 3.9 and 9.18 kyr long, were introduced at 11.5 and 34.22 m, respectively. Note that two age reversals at the top (61.9 m, 60 2 kyr) and in the middle of the section (40.28 m, 92 5 kyr) were corrected as compared to the original publications by Trauth et al. (2001, 2003) before the experiment.
monotonically increasing ages from top to bottom of the section, I have edited the age of these two layers within the errors bars of the 40 Ar/39Ar ages: 58 2 kyr and 88 5 kyr. This is quite permissible in our exercise run of the new technique, in particular since the modification of the ages does not affect the modelling result. In practice, individual ages should be eliminated as outliers rather than being edited. Table 2 lists the 40Ar/39Ar after modification of the ages at 61.9 m and 40.28 m height. The second step of the experiment is the classification of the deposits in the w60 m thick sequence in three types of deposits, (1) tephra, (2) clastic sediments, and (3) diatomite, based on the section published by Trauth et al. (2001, 2003). The sedimentation rates of these deposit types are very different, ranging from less than a tenth of a millimetre per year for diatomite to several meters of volcanic air fall deposits per within a couple of hours, followed by
Table 2 Laser-fusion 40Ar/39Ar dates and one-sigma errors from the lacustrine sequences exposed in the Ol Njorowa Gorge, south of Lake Naivasha, Kenya (Trauth et al., 2001, 2003). Ages at 61.9 m and 40.28 m height have been edited from their published values of 60 2 and 93 5 kyr within one-sigma error ranges in order to avoid age reversals. The modification of the ages does not significantly influence the result of the age modelling. Height (m) 3.43 10.20 11.65 13.13 18.00 20.92 29.22 37.45 37.79 38.09 40.28 44.86 47.48 50.28 54.83 61.90
Age (kyr)
1s error (kyr)
176.00 161.00 146.00 141.00 113.00 108.00 107.00 93.00 91.00 89.00 88.00 81.00 73.00 72.70 59.00 58.00
2.00 4.00 2.00 3.00 2.00 7.00 4.00 3.00 2.00 4.00 5.00 4.00 3.00 1.80 4.00 2.00
70
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
a longer time of no deposition. Most importantly, a thick welded tuff called the Kedong Valley Tuff, 106 4 kyr old, is located in the middle of the section, which is considered to be deposited within zero time on the 40Ar/39Ar time scale. I used a shape factor of 2 for all sediment types and a scale factor of 1.5, 2.0, and 8.0 for the three sediment types tephra, clastic sediments, and diatomite, determined in a first run of the model to adjust the Gamma parameters. Running the model with the radiometric ages and sediment types, I get the normalised Gamma/Gaussian distributions of the deposition times per sediment slice down section. From these, I can calculate the percentage overlap of the normalised distributions. Based on the percentage overlap and the information of the 10 m thick and 107 4 kyr old Kedong Valley Tuff between, possible locations of hiatuses can be determined. First, we find a significant offset of the age model and the actual radiometric ages between 10 and 11.65 m height. The percentage overlap clearly indicates a very low value of w11%, suggesting an hiatus within this height interval. An hiatus of 3.9 kyr between 151.8 and 147.8 kyr nicely compensates for the offset between the age model and the radiometric ages. Second, the zero-deposition time of the Kedong Valley Tuff is considered by setting the age of the section between 24.22 and 34.22 m to 107 kyr. Third, the correction of the age between 24.22 and 34.22 m requires the introduction of a second 9.18 kyr hiatus between 108 and 98.82 kyr at 34.22 m, i.e. on top of the welded tuff. This offset is necessary to provide a linear age model above the tuff, even if the minimum in the percentage overlap is only a local minimum. Further hiatuses may be present in the section, e.g. in the uppermost 10 m, but this is the part of the section were we have modified the age of the uppermost layer from 60 2 to 58 2 kyr, and therefore we do not touch this part of the section. Having determined the best-fitting age model with the new technique, it raises the question whether it is really better than the one created by the other methods, which were mentioned in the introduction. Undoubtedly, the statistical or numerical methods are better than visual or intuitive estimation of the age-depth relationship (Blaauw, 2010). Age models created by linear interpolation or regression between the radiometric age dates at least take into account abrupt changes in the sedimentation rates at the location of the age control points, which is suppressed by a polynomial or spline method. Indeed, both the synthetic and real example do show such abrupt changes in the sedimentation rate and therefore I consider the method better than interpolation and regression methods. The interactive or iterative adaptation of the parameters of the Gamma distribution for the unit deposition time is quite similar to the one used by the Bayesian-based BACON method by Blaauw and Christen (2011). The disadvantage of BACON in its present form, however, is the missing possibility of multiple sediment types with different Gamma distributions and the missing algorithm to detect hiatuses. The use of a memory to suppress abrupt variations in the sedimentation rates maybe suitable for radiocarbon-dated sediment cores with relatively homogeneous types of deposits. In stratigraphic sections like the one exposed in the Ol Njorowa Gorge, abrupt change in the sedimentation rates including the presence of hiatuses require a different technique than BACON, which is represented by the modelling technique presented here. 5. Conclusions Complex stratigraphic sequences with abrupt changes in the sedimentation rates, including possible hiatuses, represent major challenges for modelling the age-height relationship in paleoclimate studies. Comparing the estimates of deposition times from radiometric age determinations and unit deposition times of the sediment types included in the section helps to identify such abrupt
changes, in particular episodes of no deposition (hiatuses). The new technique is particularly useful to build an age model in the diverse fluvio-lacustrine sequences of rift lakes such the ones in Eastern Africa. Acknowledgements I thank Verena Förster and Mark Maslin for their comments on earlier versions of the manuscripts. I also thank Frank Scherbaum, Christopher Bronk Ramsey and Maarten Blaauw for inspiring discussions about Bayesian methods of age-depth modelling. I also thank the reviewers Pete Sadler and Greg Balco, as well as the editor Paul Renne for their valuable comments on the manuscript. Editorial handling by: P.R. Renne Appendix A. Supplementary data Supplementary data related to this article can be found at http:// dx.doi.org/10.1016/j.quageo.2014.03.001. References Behrensmeyer, A.K., Potts, R., Deino, A., Ditchfield, P., 2002. Olorgesailie, Kenya: a million yrs in the life of a rift basin. In: Renaut, R.W., Ashley, G.M. (Eds.), Sedimentation in Continental Rifts, SEPM Special Publication 73, pp. 97e106. Bergner, A.G.N., Trauth, M.H., 2004. Comparison of the hydrologic and hydrochemical evolution of Lake Naivasha (Kenya) during three highstands between 175 and 60 kyr BP. Palaeogeogr. Palaeoclimatol. Palaeoecol. 215 (1e2), 17e36. Blaauw, M., 2010. Methods and code for ‘classical’ age-modelling of radiocarbon sequences. Quat. Geochronol. 5, 512e518. Blaauw, M., Christen, J.A., 2011. Flexible paleoclimate age-depth models using an autoregressive gamma process. Bayesian Anal. 6, 457e474. Blaauw, M., Heegaard, E., 2012. Estimation of age-depth relationships. In: Birks, H.J.B., Juggins, S., Lotter, A., Smol, J.P. (Eds.), Data Handling and Statistical Techniques, Tracking Environmental Change Using Lake Sediments, vol. 5. Springer Verlag, Berlin, Heidelberg, New York, Tokyo, pp. 379e413. Bronk Ramsey, C., 2008. Deposition models for chronological records. Quat. Sci. Rev. 27, 42e60. Bronk Ramsey, C., 2009. Dealing with outliers and offsets in radiocarbon dating. Radicarbon 51, 1023e1045. Brown, F.H., Feibel, C.S., 1991. Stratigraphy, depositional environments, and palaeogeography of the Koobi Fora Formation. In: Harris, J.M. (Ed.), The Fossil Ungulates: Geology, Fossil Artiodactyls, and Palaeoenvironments, Koobi Fora Research Project, vol. 3. Clarendon Press, Oxford, pp. 1e30. Deino, A.L., Potts, R., 1990. Single-crystal 40Ar/39Ar dating of the Olorgesailie formation, southern Kenya rift. J. Geophys. Res. B Solid Earth Planets 95, 8453e8470. Deino, A.L., Kingston, J.D., Glen, J.M., Edgar, R.K., Hill, A., 2006. Precessional forcing of lacustrine sedimentation in the late Cenozoic Chemeron Basin, Central Kenya rift, and calibration of the Gauss/Matuyama Boundary. Earth Planet. Sci. Lett. 247, 41e60. Einsele, G., 2000. Sedimentary Basins, Evolution, Facies and Sediment Budget. Springer Verlag, Berlin, Heidelberg, New York, Tokyo. Fisher, R.V., Schmincke, H.-U., 1984. Pyroclastic Rocks. Springer Verlag, Berlin, Heidelberg, New York, Tokyo. Gasse, F., 1980. Les Diatomées lacustres Plio-Pléistocenes du Gadeb (Éthiopie): Systématique, paléoécologie, biostratigrahie. Rev. Algol. 3, 1e249. Goring, S., Williams, J.W., Blois, J.L., Jackson, S.T., Paciorek, C.J., Booth, R.K., Marlon, J.R., Blaauw, M., Christen, J.A., 2012. Deposition times in the northeastern United States during the Holocene: establishing valid priors for Bayesian age models. Quat. Sci. Rev. 48, 54e60. Hercman, H., Pawlak, J., 2012. MOD-AGE: an age-depth model construction algorithm. Quat. Geochronol. 12, 1e10. Hinderer, M., Einsele, G., 2001. The world’s large lake basins as denudationaccumulation systems and implications for their lifetimes. J. Paleolimnol. 26, 355e372. Holdship, S.A., 1976. The Paleolimnology of Lake Manyara, Tanzania: a Diatom Analysis of a 56 meter Sediment Core (Ph.D. thesis). Duke University, Durham, USA. Joordens, J.C.A., Vonhof, H.B., Feibel, C.S., Lourens, L.J., Dupont-Nivet, G., Van der Lubbe, J.H.J.L., Sier, M.J., Davies, G.R., Kroon, D., 2011. An astronomically-tuned climate framework for hominins in the Turkana Basin. Earth Planet. Sci. Lett. 307, 1e8. Kominz, M.A., Bond, G.C., 1990. A new method of testing periodicity in cyclic sediments: application to the Newark Supergroup. Earth Plenetary Sci. Lett. 98, 233e244. Kominz, M.A., Bond, G.C., 1992. Documenting the reliability and utility of the gamma method as applied to cyclic sections using forward modeling. Earth Plenetary Sci. Lett. 113, 449e457.
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71 Maher, L.J., 1972. Absolute pollen diagram of Redrock Lake, Boulder County, Colorado. Quat. Res. 2, 531e553. Owen, R.B., Crossley, R., 1992. Spatial and temporal distribution of diatoms in sediments of Lake Malawi, Central Africa, and ecological implications. J. Paleolimnol. 7, 55e71. Owen, R.B., Potts, R., Behrensmeyer, A.K., 2009. Reply to the comment on “Diatomaceous sediments and environmental change in the Pleistocene Olorgesailie Formation, southern Kenya Rift Valley” by R.B. Owen, R. Potts, A.K. Behrensmeyer and P. Ditchfield. Palaeogeogr. Palaeoclimatol. Palaeoecol. 282, 147e148. Partridge, T.C., deMenocal, P.B., Lorentz, S.A., Paiker, M.J., Vogel, J.C., 1997. Orbital forcing of climate over South Africa: a 200,000-year rainfall record from the Pretoria Saltpan. Quat. Sci. Rev. 16, 1125e1133. Sadler, P.M., 1981. Sediment accumulation rates and the completeness of stratigraphic sections. J. Geol. 89, 569e584. Sadler, P.M., 1999. The influence of hiatuses on sediment accumulation rates. GeoResearch Forum 5, 15e40. Schumer, R., Jerolmack, D.J., 2009. Real and apparent changes in sediment deposition rates through time. J. Geophys. Res. 114, F00A06. Sommerfield, C.K., 2006. On sediment accumulation rates and stratigraphic completeness: lessons from Holocene ocean margins. Cont. Shelf Res. 26, 2225e2240.
71
Steier, P., Rom, W., 2000. The use of Bayesian statistics for 14C dates of chronologically ordered samples: a critical analysis. Radiocarbon 42, 193e198. Steier, P., Rom, W., Puchegger, S., 2001. New methods and critical aspects of Bayesian mathematics for 14C calibration. Radiocarbon 43, 373e380. Trauth, M.H., Maslin, M.A., 2009. Comments on “Diatomaceous sediments and environmental change in the pleistocene Olorgesailie Formation, southern Kenya rift” by R. Bernhart Owen, Richard Potts, Anna K. Behrensmeyer and Peter Ditchfield. Palaeogeogr. Palaeoclimatol. Palaeoecol. 282, 145e146. Trauth, M.H., Deino, A., Strecker, M.R., 2001. Response of the East African climate to orbital forcing during the last interglacial (130-117 kyr BP) and the early Last Glacial (117-60 kyr BP). Geology 29, 499e502. Trauth, M.H., Deino, A., Bergner, A.G.N., Strecker, M.R., 2003. East African climate change and orbital forcing during the last 175 kyr BP. Earth Planet. Sci. Lett. 206, 297e313. Trauth, M.H., Maslin, M.A., Deino, A., Strecker, M.R., 2005. Late Cenozoic moisture history of east Africa. Science 309, 2051e2053. Trauth, M.H., Maslin, M.A., Deino, A., Strecker, M.R., Bergner, A.G.N., Dühnforth, M., 2007. High- and low-latitude forcing of Plio-Pleistocene African climate and human evolution. J. Hum. Evol. 53, 475e486.
Contents lists available at ScienceDirect
Quaternary Geochronology journal homepage: www.elsevier.com/locate/quageo
Research paper
A new probabilistic technique to build an age model for complex stratigraphic sequences Martin H. Trauth* Institut für Erd- und Umweltwissenschaften, Universität Potsdam, 14476 Potsdam, Germany
a r t i c l e i n f o
a b s t r a c t
Article history: Received 25 June 2013 Received in revised form 18 February 2014 Accepted 4 March 2014 Available online 20 March 2014
The age models of fluvio-lacustrine sedimentary sequences are often subject of discussions in paleoclimate research. The techniques employed to build an age model are very diverse, ranging from visual or intuitive estimation of the age-depth relationship over linear or spline interpolations between age control points to sophisticated Bayesian techniques also taking into account the most likely deposition times of the type of sediment within the sequence. All these methods, however, fail in detecting abrupt variations in sedimentation rates, including the possibility of episodes of no deposition (hiatus), which is the strength of the method presented in this work. The new technique simply compares the deposition time of equally thick sediment slices from the differences of subsequent radiometric age dates and the unit deposition times of the various sediment types. The percentage overlap of the distributions of these two sources of information, together with the evidence from the sedimentary record, helps to build an age model of complex sequences including abrupt variations in the rate of deposition including one or many hiatuses. Ó 2014 Elsevier B.V. All rights reserved.
Keywords: Age-depth modelling Stratigraphy Sedimentation rate Hiatus MATLAB
1. Introduction A precise age-depth model is an important prerequisite for the interpretation of high-resolution paleoclimate records. Unfortunately, these methods, although often published with downloadable computer code, guidelines and work-through examples, are hardly used in paleoclimate research (Blaauw, 2010). Many paleoclimate studies do not provide any information on the algorithm used in calculating the age model, if they ever use one, or they do not provide any information about the ambiguity and uncertainties of the result. The available age-depth modelling techniques used range from (1) visual or intuitive estimation of the age-depth relationship (e.g. Brown and Feibel, 1991; Trauth et al., 2001; Behrensmeyer et al., 2002), (2) linear, polynomial or spline interpolation or regression between the radiometric age dates (e.g. Maher, 1972; Blaauw and Heegaard, 2012), (3) calibrating the stratigraphy to insolation or orbital target curves (e.g. Partridge et al., 1997; Joordens et al., 2011), (4) forward modelling of facies variations in cyclic sections (e.g. Kominz and Bond, 1990, 1992), (5) Monte-Carlo modelling of age distributions along sediment cores (Hercman and Pawlak, 2012) to (6) Bayesian age-depth models (e.g.
* Tel.: þ49 331 977 5810; fax: þ49 331 977 5700. E-mail addresses: [email protected], [email protected]. http://dx.doi.org/10.1016/j.quageo.2014.03.001 1871-1014/Ó 2014 Elsevier B.V. All rights reserved.
Steier and Rom, 2000; Steier et al., 2001; Blaauw and Christen, 2011; Blaauw and Heegaard, 2012). The weakness of many of these techniques is mainly due to a lack of a suitable algorithm for consideration of abrupt changes in the sedimentation rate within stratigraphic sections or cores, including even the possibility of the complete lack of sedimentary layers at certain times (hiatus). Building an age model in these situations is complicated by the fact that the sedimentation rates depend strongly on the observed time resolution, with a very strong negative relationship between expected sedimentation rate and averaging time (e.g. Sadler, 1981, 1999) (Fig. 1). Such changes, however, are very common in stratigraphic sequences, for instance in the tectonically-active sedimentary basins of Eastern Africa (e.g. Brown and Feibel, 1991; Trauth et al., 2001; Behrensmeyer et al., 2002). In these basins, the sedimentation rates range from less than 0.1 m kyr1 for diatomite (a sediment composed of the skeletons of silica algae, see Table 1), and 0.1e10 m kyr1 for clastic sediments in lakes (Einsele, 2000; Hinderer and Einsele, 2001), to more than >10 m/kyr for sands (Einsele, 2000), and several meters of volcanic air fall deposits per within a couple of hours, followed by a longer time of no deposition (Fisher and Schmincke, 1984). Whereas the sedimentation rates of most sediment types are well known for lakes in Eastern Africa, the significance of hiatuses are subject passionate discussions (Trauth et al., 2005; Trauth and Maslin, 2009; Owen et al., 2009).
66
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
Fig. 1. Possible age models for three layers of clay/silt and diatomite, constrained by two radiometric age dates, 100 and 400 kyr, and their errors. The difference age models include different sedimentation rates, even with similar types of sediment, as well as one or more hiatuses of different duration. The sedimentation rates of the individual types of sediment are averaged over the entire interval. As an example, the lower clay/silt bed in Age Model 2 has a sedimentation rate of 10 m/150 kyr z 0.07 m kyr1 ¼ 0.07 m kyr1, that of the upper clay/silt bed is 10 m/50 kyr ¼ 0.2 m kyr1 ¼ 0.2 m kyr1, which averages as 20 m/200 kyr ¼ 0.1 m kyr1 ¼ 0.1 m kyr1. The choice of the age model has important implications for the sedimentation rate of the diatomite and hence for the duration of lake episode documented by the diatomite.
The possibility of an abrupt change in the sedimentation rate was even regarded as a weakness of a linear over a spline age model, even though these changes exist but are often excluded or ignored in stratigraphic sections, similar to the possibility of hiatuses. Instead, an age model has to be smooth, often by introducing an arbitrary chosen memory in the sedimentation rates along the core, avoiding extreme variations or extremely low or high sedimentation rates, without abrupt changes in the measures of central tendency and dispersion, and without major gaps (e.g. Bronk Ramsey, 2008, 2009; Blaauw and Christen, 2011). Furthermore, it is assumed that a single sediment type, defined by the Gammadistributed sedimentation rate, occurs through the section, which could be a good model in some lacustrine settings, but is an undue assumption in many sedimentary environments (e.g. Bronk Ramsey, 2008, 2009; Blaauw and Christen, 2011). Here, I present a new technique to determine the age-depth relationship for stratigraphic sections and sediment cores with
extremely fluctuating sedimentation rates, including hiatuses. The technique uses two independent sources of information to build an age model, (1) the time of deposition based on radiometric ages with their Gaussian errors and (2) the typical unit deposition time of the various sediment types (Fig. 2). Although the algorithm has been designed, tested and applied to age dates with Gaussian errors, it can be adapted for other applications quite easily. The agreement, or disagreement, between the two estimates of the time of deposition from two statistically independent sources of information helps to build an age-depth relationship in complex stratigraphic sequences. While this new technique in the case of high-resolution radiocarbon chronologies of relatively young (<104 yrs) and homogeneous sediments is not as straight forward as Bayesian techniques, it is superior in older 40Ar/39Ar dated sedimentary sequences with extremely fluctuating sedimentation rates and high probabilities for the occurrence of hiatuses. 2. The system to be modelled
Table 1 Typical sedimentation rates of diatomite as documented in Quaternary lake sediment sequences from East Africa. Location Gadeb, Ethiopia Ol Njorowa Gorge Ol Njorowa Gorge Gicheru Olorgesailie Gicheru Lake Malawi Kariandusi Lake Manyara Barsemoi
Sedimentation rate 1
0.1 m kyr 0.03e0.25 m kyr1 0.13e0.4 m kyr1 0.5 m kyr1 1.0 m kyr1 1.2 m kyr1 >1.3 m kyr1 1.5 m kyr1 1.9 m kyr1 2.2 (0.3e12) m kyr1
Reference Gasse 1980 Bergner and Trauth, 2004 Trauth et al., 2001 Trauth et al., 2007 Deino and Potts, 1990 Trauth et al., 2007 Owen and Crossley, 1992 Trauth et al., 2007 Holdship, 1976 Deino et al., 2006
Sedimentary sections consist of a layered sequence of materials with different grain size and composition. The rate of deposition of these materials is described by the sedimentation rate (sediment thickness per unit time) or, alternatively, the accumulation rate (solid sediment mass per unit area and time, in kg m2 kyr1) (Einsele, 2000). The sedimentation rate (in m kyr1) is determined from the sedimentary unit of the thickness z, divided by the time difference t between two radiometric age dates below and above the unit. While the terms sedimentation and accumulation rate are usually used correctly in the literature, there is no consensus on the name of the inverse of the sedimentation rate (time per sediment thickness). I prefer the term unit deposition time (in kyr m1), which
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
Fig. 2. (A) Gaussian frequency distribution of the time difference per slice between two age dates with Gaussian errors is determine by radiometric dating. (B) The Gamma distributions of the unit deposition time defined for each of the type of deposit (e.g. tephra, clastic sediments, and diatomite) is determined by the sediment record.
is very close to deposition time used by Goring et al. (2012) but different from accumulation rate used by Blaauw and Christen (2011). The sedimentation rate varies dramatically between different depositional environments, ranging from >104 m kyr1 in continental slopes, alluvial fans and delta plains to <101 m kyr1 in lakes and the deep sea below the carbonate compensation depth (Einsele, 2000). Of course, these rates also vary within the same environment, even at the same place, driven by volcano-tectonic, climate and other factors, especially in tectonically active regions (Einsele, 2000). Unfortunately, the stratigraphic record is fragmentary, i.e. it contains stratigraphic gaps without deposition, often without any information on the exact duration of these episodes of no deposition or hiatuses (Sadler, 1981, 1999; Sommerfield, 2006; Schumer and Jerolmack, 2009). The possibility of fluctuating sedimentation rates and several hiatuses between two radiometric age dates complicates the assessment of the actual time of deposition of a sedimentary sequence. As an example, let us assume a simple sequence of three sediment layers, a unit of clastic sediments overlain by a diatomite bed, which is again overlain by clastic sediments, each of these sediment layers being 10 m thick (Fig. 1). The sedimentary sequence is bounded by two age dates, 100 and 400 kyr old. Fig. 1 shows seven possible age models for the sedimentary sequence, all of them consistent with the plausible sedimentation rates of fine-grained sediments (0.0697e0.1 m kyr1) and diatomite (0.1e10 m kyr1). Age models 3, 4 and 6, however, also have hiatuses with a duration of 40e140 kyr, which complicates the interpretation of the stratigraphic sequence. The hiatuses can occur both at the layer boundaries (age model 6) as well as within layers (age models 3 and 4). How can we find the best (or best-fitting) age model for such sedimentary sequences? In statistics, or in sciences in general, the golden rule is always to use the simplest model, which is explains the data. In mathematical terms, the model with the smallest
67
number of adjustable parameters is considered the best (or bestfitting) as long as the mean least-squared error is below a critical value. As an example, if a simple model with a linear sedimentation rate (two parameters, slope and y-intercept) describes the data well in a least-squared sense, we prefer it over a model with three parameters (cubic polynomial, three parameters). In practice, a piecewise linear or cubic spline interpolation between age control points is used to find the best age model (e.g. Blaauw and Heegaard, 2012). In this context, it has often been argued that cubic spline interpolation between age control points is better than simple linear interpolation because it smoothes out minor errors in age data points (Maher, 1972; Blaauw and Heegaard, 2012). This implies that abrupt changes in sedimentation rates are due to errors in the age dates. On the other hand, in more complicated sequences with varying sediment types, such abrupt changes in the sedimentation rates may actually exist, as shown in Fig. 1, also showing zero sedimentation rate (hiatuses); even negative sedimentation rates (erosion) are possible. In fact, the advantage of the technique presented here is the ability to detect such episodes without deposition in the record without any prior assumptions. Furthermore it allows to include multiple sediment types with different sedimentation rates. Most importantly, there are no a priori assumptions made regarding a memory or coherence in the sedimentation rates in the age model to suppress abruptly changing rates; abruptly changing rates are explicitly expected at discontinuities between layers of different sediment types. In such an age model, a linear interpolation between the age control points could provide the better age model over cubic splines. 3. Description of the algorithm I use a new technique to build an age model for stratigraphic sections and sediment cores (Fig. 1). The basic idea of the new method is to estimate the most likely time of deposition between two radiometric age dates, such as 40Ar/39Ar dates, based on two sources of information (Fig. 2). First, the time difference is determined from the two age dates and their independent Gaussian errors. Second, the time difference is again calculated from the time accumulation rate of the sediment types and their Gammadistributed dispersion. The time difference per slice between two age dates is determined by radiometric dating of volcanic ash layers in the case of 40Ar/39Ar dates of feldspar crystals or any other dating technique (Fig. 2a). The age dates are statistically independent from the type of sediment (or any other deposit) between the dated levels in the stratigraphic section. The errors of the ages are assumed to be independent from each other. Therefore, the time difference between two age dates can be calculated, for instance by using the Gaussian law of error propagation. In practice, calculating the difference of two ages in a stratigraphic section is slightly more difficult than described here. Radiometric age dates such as 40 Ar/39Ar ages of subsequent sediment layers, or more precisely, the one of the intercalated ash layers, are often not significantly different within the error bars (zero age difference) or the estimated age of the stratigraphically higher layer is older than the lower one (age reversal). For simplicity in our experiment, we assume monotonically increasing ages with depth and therefore exclude those ages that are unexpectedly younger or older than expected from the stratigraphic sequence. Of course we have to explain these outliers or design a more complex age model for the section, i.e. considering reworking of certain units in the section. The unit deposition time is defined for each of the type of deposit in the stratigraphic section (Fig. 2b). First, the type of sediment is defined in the outcrop or in the sediment core. We can work with
68
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
any number of sediment types, but at the beginning we should work with a small number. Second, a Gamma model is used for the most likely unit deposition time for each sediment type (e.g. three types of materials: tephra, clastic sediments, diatomite). The Gamma distribution has two parameters, shape and scale. In the experiments described here, a Gamma model gamma(shape, scale) with shape ¼ 2 best describes the frequency distribution of unit deposition times, which is in agreement with the work of Blaauw and Christen (2011). At the lower end of unit deposition times is diatomite with typical sedimentation rates between 0.03 and 12 m kyr1, with a mode at around 1 m kyr1 (Table 1). In our experiments, a Gamma model with gamma(2,1.5) best describes the expected range of unit deposition times. For clastic sediments we found gamma(2,2) as being the best-fitting model, corresponding to the most frequent sedimentation rate of 0.5 m kyr1. The best Gamma model for tephra is gamma(2,8) with a very wide range of possible unit deposition times, taking into account the large variation in possible sedimentation rate between several millimetres of fine-grained ash and thick layers of air fall tephra deposited within hours. The agreement, or disagreement, between the two estimates of the time of deposition from two statistically independent sources of information helps to build an age-depth relationship in complex stratigraphic sequences. If the two estimates of the time of deposition match, the probability of a hiatuses is very low. The time difference between the two estimates is regarded as the duration of the hiatuses. Because both values are linked to their frequency distributions, we can also calculate the error probability for the conclusion “hiatus” and its duration. While this new technique in the case of high-resolution radiocarbon chronologies of relatively young and homogeneous sediments is not as straight forward as Bayesian techniques, it is superior in sedimentary sequences with fluctuating sedimentation rates and high probabilities for the occurrence of hiatuses. Alternatively, if there is no evidence for an hiatus between the age dates, we can multiply the two independent frequency distributions, Gaussian for the age difference and Gamma for the sediment types, to calculate a better estimate of the time of deposition and its frequency distribution. After having determined the distribution of the deposition times along the stratigraphic record, the paleoclimatologist wants to know the best-fitting age model at which one can interpolate the physical, chemical and biological environmental indicators. Much has been written about using the mean, weighted mean, median, mode, or any other more complex measure of central tendency to build the age model from the frequency distribution. I found the argument of Blaauw and Heegaard (2012) very well justified, recommending to use the mode (modal value) of the probability distribution over any other measure of central tendency. According to them, this has the advantage of being the age that is most likely, but it can concentrate emphasis on one peak that happens to be fractionally higher than another in a distribution that is bimodal or even multi-modal, which could certainly be the case in radiocarbon-dated sections rather than in Argon-dated sequences. Similar to them, I use percentiles of the probability distribution to describe the spread or uncertainty used as tolerance, for example, 25% below and above the mode in the example below. 4. Examples 4.1. Synthetic Argon-dated lake-sediment sequence Here, I demonstrate the use of the new technique to build an age model for a synthetic 50 m thick lake sediment sequence (Fig. 3). The Supplementary material of the paper contains a MATLAB script called script.m and the input files input_section.txt and
Fig. 3. Best-fitting age model of a synthetic section with three sediment types (tephra, clastic sediments, and diatomite) and radiometric age dates (in black, with one-sigma error bars). The age model (blue line) was calculated from the modes of the multiplied distributions of the age differences per slice and the unit deposition times of the deposits, displayed together with the 25% confidence bounds below and above the modes (red lines). Based on the initial offset of the age model and the radiometric age dates above 37.4 m height and the relatively low percentage overlap (w7.4%) of the distributions, a 2 kyr long hiatus was introduced at 35 m. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
input_ages.txt as an example to the application of the technique. The MATLAB script contains all necessary steps to build the age model based on the lithological classification and ages of the sedimentary sequence in the input files. The input file input_ section.txt has the format
0 6 12 ð.Þ 38 42 47
6 12 14
1 3 1
42 47 50
1 2 1
where the first and second column include the base and top of each sedimentary unit in meters. The third column includes an identifier for the type of deposit, in which 1 stands for tephra, 2 for clastic sediments and 3 for diatomite. The input file input_ages.txt contains the Argon ages in the format
2:8 398:8 2:0 7:6 394:0 1:5 11 387:4 2:0 ð.Þ 37:4 356:7 3:0 43:2 352:8 1:6 47:8 348:8 1:2 where the first column includes the height of the sample in the section in meters. The second and third column contains the ages before present and their one-sigma errors. After loading the data into the workspace of MATLAB, we define a starting value for the maximum age maxage of the section. The initial choice of the value of maxage is not critical to the result of the experiment and can be modified in a later run of the model to get better results. The section is then divided into equally thick slices of sediment, for which the deposition time is determined from the radiometric ages and the unit deposition time of the sediment types. The choice of the thickness slices defines the observed time resolution and
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
therefore the amplitude of possible variations in the unit deposition times. Furthermore, smaller values for the thickness slices produce more accurate results at the expense of higher computation time. In the example, I used slices ¼ 0.1 for 0.1 m thick slices to get a good result. The following set of commands in the MATLAB script includes the interpolation of the sedimentary data from input_section.txt upon an evenly spaced grid of the 0.1 m slices. The first estimate of the deposition time is based on the unit deposition time of the sediment types in the section. We first design the Gamma models for the three sediment types in the section as classified in input_section.txt. Similar to the Bayesian techniques by Blaauw and Christen (2011), we can start with an initial set of shape and scale factors and then we iteratively improve the factors to get better results. In the synthetic experiment, I used a shape factor of 2 for all sediment types and a scale factor of 0.5, 0.8, and 2.0 for the three sediment types tephra, clastic sediments, and diatomite, respectively, as defined in the array shape ¼ [2 0.5;2 0.8;2 2]. The Gamma models for the sediment types helps to calculate a series of unit deposition times along the section, from which we calculate the modes and the quartiles, as well as the cumulative deposition time of the section. The second estimate of the deposition time comes from the difference of the Argon ages in the file input_ages.txt. After calculating the difference of subsequent ages, we interpolate a series of deposition ages along the section for each sediment slice, calculate the modes and the quartiles, as well as the cumulative deposition time of the section. Since we now have two independent estimates for the deposition time including the normalised Gamma/Gaussian distributions of these estimates, we can calculate the percentage overlap of the normalised distributions. The value of 1eoverlap helps to define intervals, which probably contain one or several hiatuses. If there is no reason to prefer one of the two estimates for the deposition time, then we can multiply the two distributions to obtain an improved estimate of the true deposition time for each sediment slice. In our example, we notice very low values (w7.4%) of the percentage overlap of the two distributions between 33.9 and 36 m. Here, the deposition time from the sediment types is 0.33, whereas the deposition time from the radiometric ages is 0.37 years. In contrast, the deposition time at around 10 m height is in good agreement, both being w0.19 years per 0.1 m thick slice. To compensate for the discrepancy between the two values for the deposition time, we can now install a hiatus between 33.9 and 36 m in the section, up to the full difference between the ages and the expected unit deposition time of the sediment. The exact location of the hiatus between 33.9 and 36 m is unknown, but may be found by looking into the stratigraphic logs of the outcrop. In this synthetic example, I have decided to introduce a 2 kyr long hiatus in the middle of the interval at 35 m, which gives a good age-depth model for the section. 4.2. Argon-dated lake-sediment sequence of Lake Naivasha Here, I use the new technique to build an age model for the lake sediment sequences, ca. 150e60 kyr old, exposed in the Ol Njorowa Gorge south of Lake Naivasha in the Central Kenya Rift (Trauth et al., 2001, 2003) (Fig. 4). Anorthoclase and sanidine phenocryst concentrates from 16 tephra beds were dated by the laser-fusion 40 Ar/39Ar method. The chronology is free of age reversals within the one-sigma error bars of the 40Ar/39Ar ages. The means, however, show in fact two reversals at the top (61.9 m, 60 2 kyr) and in the middle of the section (40.28 m, 92 5 kyr). Since my technique, similar to BACON by Blaauw and Christen (2011), requires
69
Fig. 4. Best-fitting age model of the Late Pleistocene lacustrine section exposed in the Ol Njorowa Gorge, south of Lake Naivasha, Kenya. The w60 m thick sequence consists of three types of deposits, (1) tephra, (2) clastic sediments, and (3) diatomite, based on the stratigraphic column of the section published by Trauth et al. (2001, 2003). Based on the percentage overlap of the distributions of the age differences per slice and the unit deposition times of the deposits, also considering the presence of the 10 m thick Kedong Valley Tuff (107 kyr old) in the middle of the section, two hiatuses, 3.9 and 9.18 kyr long, were introduced at 11.5 and 34.22 m, respectively. Note that two age reversals at the top (61.9 m, 60 2 kyr) and in the middle of the section (40.28 m, 92 5 kyr) were corrected as compared to the original publications by Trauth et al. (2001, 2003) before the experiment.
monotonically increasing ages from top to bottom of the section, I have edited the age of these two layers within the errors bars of the 40 Ar/39Ar ages: 58 2 kyr and 88 5 kyr. This is quite permissible in our exercise run of the new technique, in particular since the modification of the ages does not affect the modelling result. In practice, individual ages should be eliminated as outliers rather than being edited. Table 2 lists the 40Ar/39Ar after modification of the ages at 61.9 m and 40.28 m height. The second step of the experiment is the classification of the deposits in the w60 m thick sequence in three types of deposits, (1) tephra, (2) clastic sediments, and (3) diatomite, based on the section published by Trauth et al. (2001, 2003). The sedimentation rates of these deposit types are very different, ranging from less than a tenth of a millimetre per year for diatomite to several meters of volcanic air fall deposits per within a couple of hours, followed by
Table 2 Laser-fusion 40Ar/39Ar dates and one-sigma errors from the lacustrine sequences exposed in the Ol Njorowa Gorge, south of Lake Naivasha, Kenya (Trauth et al., 2001, 2003). Ages at 61.9 m and 40.28 m height have been edited from their published values of 60 2 and 93 5 kyr within one-sigma error ranges in order to avoid age reversals. The modification of the ages does not significantly influence the result of the age modelling. Height (m) 3.43 10.20 11.65 13.13 18.00 20.92 29.22 37.45 37.79 38.09 40.28 44.86 47.48 50.28 54.83 61.90
Age (kyr)
1s error (kyr)
176.00 161.00 146.00 141.00 113.00 108.00 107.00 93.00 91.00 89.00 88.00 81.00 73.00 72.70 59.00 58.00
2.00 4.00 2.00 3.00 2.00 7.00 4.00 3.00 2.00 4.00 5.00 4.00 3.00 1.80 4.00 2.00
70
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71
a longer time of no deposition. Most importantly, a thick welded tuff called the Kedong Valley Tuff, 106 4 kyr old, is located in the middle of the section, which is considered to be deposited within zero time on the 40Ar/39Ar time scale. I used a shape factor of 2 for all sediment types and a scale factor of 1.5, 2.0, and 8.0 for the three sediment types tephra, clastic sediments, and diatomite, determined in a first run of the model to adjust the Gamma parameters. Running the model with the radiometric ages and sediment types, I get the normalised Gamma/Gaussian distributions of the deposition times per sediment slice down section. From these, I can calculate the percentage overlap of the normalised distributions. Based on the percentage overlap and the information of the 10 m thick and 107 4 kyr old Kedong Valley Tuff between, possible locations of hiatuses can be determined. First, we find a significant offset of the age model and the actual radiometric ages between 10 and 11.65 m height. The percentage overlap clearly indicates a very low value of w11%, suggesting an hiatus within this height interval. An hiatus of 3.9 kyr between 151.8 and 147.8 kyr nicely compensates for the offset between the age model and the radiometric ages. Second, the zero-deposition time of the Kedong Valley Tuff is considered by setting the age of the section between 24.22 and 34.22 m to 107 kyr. Third, the correction of the age between 24.22 and 34.22 m requires the introduction of a second 9.18 kyr hiatus between 108 and 98.82 kyr at 34.22 m, i.e. on top of the welded tuff. This offset is necessary to provide a linear age model above the tuff, even if the minimum in the percentage overlap is only a local minimum. Further hiatuses may be present in the section, e.g. in the uppermost 10 m, but this is the part of the section were we have modified the age of the uppermost layer from 60 2 to 58 2 kyr, and therefore we do not touch this part of the section. Having determined the best-fitting age model with the new technique, it raises the question whether it is really better than the one created by the other methods, which were mentioned in the introduction. Undoubtedly, the statistical or numerical methods are better than visual or intuitive estimation of the age-depth relationship (Blaauw, 2010). Age models created by linear interpolation or regression between the radiometric age dates at least take into account abrupt changes in the sedimentation rates at the location of the age control points, which is suppressed by a polynomial or spline method. Indeed, both the synthetic and real example do show such abrupt changes in the sedimentation rate and therefore I consider the method better than interpolation and regression methods. The interactive or iterative adaptation of the parameters of the Gamma distribution for the unit deposition time is quite similar to the one used by the Bayesian-based BACON method by Blaauw and Christen (2011). The disadvantage of BACON in its present form, however, is the missing possibility of multiple sediment types with different Gamma distributions and the missing algorithm to detect hiatuses. The use of a memory to suppress abrupt variations in the sedimentation rates maybe suitable for radiocarbon-dated sediment cores with relatively homogeneous types of deposits. In stratigraphic sections like the one exposed in the Ol Njorowa Gorge, abrupt change in the sedimentation rates including the presence of hiatuses require a different technique than BACON, which is represented by the modelling technique presented here. 5. Conclusions Complex stratigraphic sequences with abrupt changes in the sedimentation rates, including possible hiatuses, represent major challenges for modelling the age-height relationship in paleoclimate studies. Comparing the estimates of deposition times from radiometric age determinations and unit deposition times of the sediment types included in the section helps to identify such abrupt
changes, in particular episodes of no deposition (hiatuses). The new technique is particularly useful to build an age model in the diverse fluvio-lacustrine sequences of rift lakes such the ones in Eastern Africa. Acknowledgements I thank Verena Förster and Mark Maslin for their comments on earlier versions of the manuscripts. I also thank Frank Scherbaum, Christopher Bronk Ramsey and Maarten Blaauw for inspiring discussions about Bayesian methods of age-depth modelling. I also thank the reviewers Pete Sadler and Greg Balco, as well as the editor Paul Renne for their valuable comments on the manuscript. Editorial handling by: P.R. Renne Appendix A. Supplementary data Supplementary data related to this article can be found at http:// dx.doi.org/10.1016/j.quageo.2014.03.001. References Behrensmeyer, A.K., Potts, R., Deino, A., Ditchfield, P., 2002. Olorgesailie, Kenya: a million yrs in the life of a rift basin. In: Renaut, R.W., Ashley, G.M. (Eds.), Sedimentation in Continental Rifts, SEPM Special Publication 73, pp. 97e106. Bergner, A.G.N., Trauth, M.H., 2004. Comparison of the hydrologic and hydrochemical evolution of Lake Naivasha (Kenya) during three highstands between 175 and 60 kyr BP. Palaeogeogr. Palaeoclimatol. Palaeoecol. 215 (1e2), 17e36. Blaauw, M., 2010. Methods and code for ‘classical’ age-modelling of radiocarbon sequences. Quat. Geochronol. 5, 512e518. Blaauw, M., Christen, J.A., 2011. Flexible paleoclimate age-depth models using an autoregressive gamma process. Bayesian Anal. 6, 457e474. Blaauw, M., Heegaard, E., 2012. Estimation of age-depth relationships. In: Birks, H.J.B., Juggins, S., Lotter, A., Smol, J.P. (Eds.), Data Handling and Statistical Techniques, Tracking Environmental Change Using Lake Sediments, vol. 5. Springer Verlag, Berlin, Heidelberg, New York, Tokyo, pp. 379e413. Bronk Ramsey, C., 2008. Deposition models for chronological records. Quat. Sci. Rev. 27, 42e60. Bronk Ramsey, C., 2009. Dealing with outliers and offsets in radiocarbon dating. Radicarbon 51, 1023e1045. Brown, F.H., Feibel, C.S., 1991. Stratigraphy, depositional environments, and palaeogeography of the Koobi Fora Formation. In: Harris, J.M. (Ed.), The Fossil Ungulates: Geology, Fossil Artiodactyls, and Palaeoenvironments, Koobi Fora Research Project, vol. 3. Clarendon Press, Oxford, pp. 1e30. Deino, A.L., Potts, R., 1990. Single-crystal 40Ar/39Ar dating of the Olorgesailie formation, southern Kenya rift. J. Geophys. Res. B Solid Earth Planets 95, 8453e8470. Deino, A.L., Kingston, J.D., Glen, J.M., Edgar, R.K., Hill, A., 2006. Precessional forcing of lacustrine sedimentation in the late Cenozoic Chemeron Basin, Central Kenya rift, and calibration of the Gauss/Matuyama Boundary. Earth Planet. Sci. Lett. 247, 41e60. Einsele, G., 2000. Sedimentary Basins, Evolution, Facies and Sediment Budget. Springer Verlag, Berlin, Heidelberg, New York, Tokyo. Fisher, R.V., Schmincke, H.-U., 1984. Pyroclastic Rocks. Springer Verlag, Berlin, Heidelberg, New York, Tokyo. Gasse, F., 1980. Les Diatomées lacustres Plio-Pléistocenes du Gadeb (Éthiopie): Systématique, paléoécologie, biostratigrahie. Rev. Algol. 3, 1e249. Goring, S., Williams, J.W., Blois, J.L., Jackson, S.T., Paciorek, C.J., Booth, R.K., Marlon, J.R., Blaauw, M., Christen, J.A., 2012. Deposition times in the northeastern United States during the Holocene: establishing valid priors for Bayesian age models. Quat. Sci. Rev. 48, 54e60. Hercman, H., Pawlak, J., 2012. MOD-AGE: an age-depth model construction algorithm. Quat. Geochronol. 12, 1e10. Hinderer, M., Einsele, G., 2001. The world’s large lake basins as denudationaccumulation systems and implications for their lifetimes. J. Paleolimnol. 26, 355e372. Holdship, S.A., 1976. The Paleolimnology of Lake Manyara, Tanzania: a Diatom Analysis of a 56 meter Sediment Core (Ph.D. thesis). Duke University, Durham, USA. Joordens, J.C.A., Vonhof, H.B., Feibel, C.S., Lourens, L.J., Dupont-Nivet, G., Van der Lubbe, J.H.J.L., Sier, M.J., Davies, G.R., Kroon, D., 2011. An astronomically-tuned climate framework for hominins in the Turkana Basin. Earth Planet. Sci. Lett. 307, 1e8. Kominz, M.A., Bond, G.C., 1990. A new method of testing periodicity in cyclic sediments: application to the Newark Supergroup. Earth Plenetary Sci. Lett. 98, 233e244. Kominz, M.A., Bond, G.C., 1992. Documenting the reliability and utility of the gamma method as applied to cyclic sections using forward modeling. Earth Plenetary Sci. Lett. 113, 449e457.
M.H. Trauth / Quaternary Geochronology 22 (2014) 65e71 Maher, L.J., 1972. Absolute pollen diagram of Redrock Lake, Boulder County, Colorado. Quat. Res. 2, 531e553. Owen, R.B., Crossley, R., 1992. Spatial and temporal distribution of diatoms in sediments of Lake Malawi, Central Africa, and ecological implications. J. Paleolimnol. 7, 55e71. Owen, R.B., Potts, R., Behrensmeyer, A.K., 2009. Reply to the comment on “Diatomaceous sediments and environmental change in the Pleistocene Olorgesailie Formation, southern Kenya Rift Valley” by R.B. Owen, R. Potts, A.K. Behrensmeyer and P. Ditchfield. Palaeogeogr. Palaeoclimatol. Palaeoecol. 282, 147e148. Partridge, T.C., deMenocal, P.B., Lorentz, S.A., Paiker, M.J., Vogel, J.C., 1997. Orbital forcing of climate over South Africa: a 200,000-year rainfall record from the Pretoria Saltpan. Quat. Sci. Rev. 16, 1125e1133. Sadler, P.M., 1981. Sediment accumulation rates and the completeness of stratigraphic sections. J. Geol. 89, 569e584. Sadler, P.M., 1999. The influence of hiatuses on sediment accumulation rates. GeoResearch Forum 5, 15e40. Schumer, R., Jerolmack, D.J., 2009. Real and apparent changes in sediment deposition rates through time. J. Geophys. Res. 114, F00A06. Sommerfield, C.K., 2006. On sediment accumulation rates and stratigraphic completeness: lessons from Holocene ocean margins. Cont. Shelf Res. 26, 2225e2240.
71
Steier, P., Rom, W., 2000. The use of Bayesian statistics for 14C dates of chronologically ordered samples: a critical analysis. Radiocarbon 42, 193e198. Steier, P., Rom, W., Puchegger, S., 2001. New methods and critical aspects of Bayesian mathematics for 14C calibration. Radiocarbon 43, 373e380. Trauth, M.H., Maslin, M.A., 2009. Comments on “Diatomaceous sediments and environmental change in the pleistocene Olorgesailie Formation, southern Kenya rift” by R. Bernhart Owen, Richard Potts, Anna K. Behrensmeyer and Peter Ditchfield. Palaeogeogr. Palaeoclimatol. Palaeoecol. 282, 145e146. Trauth, M.H., Deino, A., Strecker, M.R., 2001. Response of the East African climate to orbital forcing during the last interglacial (130-117 kyr BP) and the early Last Glacial (117-60 kyr BP). Geology 29, 499e502. Trauth, M.H., Deino, A., Bergner, A.G.N., Strecker, M.R., 2003. East African climate change and orbital forcing during the last 175 kyr BP. Earth Planet. Sci. Lett. 206, 297e313. Trauth, M.H., Maslin, M.A., Deino, A., Strecker, M.R., 2005. Late Cenozoic moisture history of east Africa. Science 309, 2051e2053. Trauth, M.H., Maslin, M.A., Deino, A., Strecker, M.R., Bergner, A.G.N., Dühnforth, M., 2007. High- and low-latitude forcing of Plio-Pleistocene African climate and human evolution. J. Hum. Evol. 53, 475e486.