Interactive comment on “ Interpolation methods for Antarctic ice-core timescales : application to Byrd , Siple Dome and Law Dome ice cores ” by T . J .

Abstract. Antarctic ice cores have often been dated by matching distinctive features of atmospheric methane to those detected in annually dated ice cores from Greenland. Establishing the timescale between these tie-point ages requires interpolation. While the uncertainty at tie points is relatively well described, uncertainty of the interpolation is not. Here we assess the accuracy of three interpolation schemes using data from the WAIS Divide ice core in West Antarctica; we compare the interpolation methods with the annually resolved timescale for the past 30 kyr. Linear interpolation yields large age errors (up to 380 years) between tie points, abrupt changes in duration of climate events at tie points, and an age bias. Interpolations based on the smoothest accumulation rate (ACCUM) or the smoothest annual-layer thickness (ALT) yield timescales that more closely agree with the annually resolved timescale and do not have abrupt changes in duration at tie points. We use ALT to assess the uncertainty in existing timescales for the past 30 kyr from Byrd, Siple Dome, and Law Dome. These ice-core timescales were developed with methods similar to linear interpolation. Maximum age differences exceed 1000 years for Byrd and Siple Dome, and 500 years for Law Dome. For the glacial–interglacial transition (21 to 12 kyr), the existing timescales are, on average, older than ALT by 40 years for Byrd, 240 years for Siple Dome, and 150 years for Law Dome. Because interpolation uncertainty is often not considered, age uncertainties for ice-core records are often underestimated.


Introduction
Interpretation of paleoclimate records depends on accurate chronologies.Ice cores provide exceptional records of past climate (Grootes et al., 1993;EPICA Members, 2006) and are among the best-dated paleoclimate records (Meese et al., 1997;Svensson et al., 2008).In Greenland, the icecore timescales are based on identification of annual layers preserved in the ice; the Greenland Ice Core Chronology (GICC05, Vinther et al., 2006;Rasmussen et al., 2006;Andersen et al., 2006;Svensson et al., 2008) has identified layers to 60 kyr (thousands of years before 1950).In Antarctica, dating of ice cores has been challenging because low accumulation at many sites hampers annual-layer identification.Thus, Antarctic timescales are often derived by matching distinctive age markers detected and dated in Greenland timescales.Abrupt variations of atmospheric methane are the most commonly used tie points (Blunier et al., 1998), and some combination of approximately 11 events in the past 30 kyr (Fig. 1) have been used to transfer Greenland timescales to Antarctic ice cores.Other potential markers such as volcanic ash, cosmogenic isotope anomalies (Raisbeck et al., 2007), and sulfate peaks (Sigl et al., 2013;Svensson et al., 2013) either are more difficult to identify, are less spatially extensive, or occur less frequently.The methane tie points allow the gases to be dated, but there is an additional step to derive the ice timescale.For any given depth in an ice core, the gas trapped in bubbles is younger than the ice because the gas is not trapped until 50 to 100 m below the surface.The difference between the age of the ice and the age of the gas, termed age, must be added to gas ages to derive the ice timescale.
Published by Copernicus Publications on behalf of the European Geosciences Union.The age uncertainty at the methane tie points is well described (e.g., Blunier et al., 1998Blunier et al., , 2007;;Pedro et al., 2011;Stenni et al., 2011), but uncertainty introduced by the interpolation scheme is often given less attention.Two common approaches have been used to interpolate between methane tie points.In the first approach, variations of linear interpolation were used to construct timescales for Taylor Dome (Steig et al., 1998), Byrd (Blunier and Brook, 2001), and Siple Dome (Brook et al., 2005).Law Dome (Pedro et al., 2011) does not use linear interpolation, instead assuming constant accumulation rates between tie points and a thinning function computed with a one-dimensional ice-flow model.The resulting timescale, however, has many of the same features as a timescale constructed with linear interpolation, as is shown below.In the second approach, Bayesian statistical methods are applied to multiple ice cores simultaneously (Lemieux-Dudon et al., 2010).This method starts with initial timescales derived from ice-flow modeling driven by accumulation-rate histories based on the stable-isotope records.The initial timescales are then adjusted to optimize the agreement among initial timescales and age constraints by taking into account their respective confidence intervals.This method has been used to construct consistent timescales for EDML, EDC, Vostok, TALDICE, and NGRIP (Lemieux-Dudon et al., 2010;Buiron et al., 2011;Schüpbach et al., 2011;Veres et al., 2013;Bazin et al., 2013).Uncertainties are computed based on the uncertainties of tie points and variances of the accumulation rate, thinning function, and lock-in depth used in the initial timescale.
Assessing the interpolation uncertainty of previous Antarctic timescales is not straightforward because there is no "true" timescale with which to compare.The new WAIS Divide timescale, with annual resolution for the past 30 kyr, provides an opportunity to assess the interpolation uncertainty; by constructing timescales for WAIS Divide assuming known tie-point ages, the various interpolation methods can be compared to the "true" annually resolved timescale.We focus on three interpolation methods: linear, smoothest accumulation rate (ACCUM), and smoothest annual-layer thickness (ALT).ACCUM is based on the work of Waddington et al. (2003), which showed that linear interpolations of depthage between known tie points imply artificial saw-tooth histories in accumulation rate, even in the presence of steadystate ice flow.In order to get a physically based depthage interpolation, they introduced an inverse procedure that inferred a smooth accumulation-rate history such that the depth-age calculated with an ice-flow model matched the tie points to an acceptable tolerance (e.g., Lundin, 2012).ALT was developed as a complement to ACCUM for situations where the layer thinning is not well constrained by ice-flow modeling.It was inspired by annual-layer interpretations when the best choice for an uncertain year is the one that matches the average annual-layer thickness best.
In this work, we use ACCUM and ALT in the limiting case where tie-point ages are assumed to be known exactly.This allows the interpolation uncertainty to be isolated from the tie-point uncertainty.We show that ACCUM and ALT agree with the WAIS Divide annual timescale better than linear interpolation, and then we discuss the uncertainty of three existing timescales (Siple Dome, Byrd, and Law Dome) that use near-linear interpolation.

Methods
We evaluate interpolation methods using the WDC06A-7 timescale for the WAIS Divide ice core (WAIS Divide Project Members, 2013).The WDC06A-7 timescale is based primarily on electrical measurements and is annually resolved to 29.6 kyr (2800 m depth of a total ice to thickness of ∼ 3450 m).We use WDC06A-7 as the "true" timescale.Eleven ice-age tie points (Table 1) are used to evaluate three interpolation schemes by interpolating between these tie points and comparing results with the "true" timescale.We use tie points similar to what has been used in other Antarctic timescales (e.g., ages of identifiable variations in methane) except during the mid-and late Holocene when they are evenly spaced (2, 4, and 6 kyr).An additional tie point older than the annual timescale (35.8 kyr) is included to constrain the variations in annual-layer thickness at the older boundary.
The interpolation methods are 1.Linear -yields constant annual-layer thicknesses between age markers 2. ACCUM -uses a simple ice-flow model and an inverse method to find the smoothest accumulation-rate history that matches the depth-age markers 3. ALT -uses an inverse method to find the smoothest progression of annual-layer thicknesses that fit the depthage markers.
In all cases, the interpolations are forced to match the age of tie points nearly exactly.We do not use the stable-isotope record as a guide for interpolation.At many Antarctic ice-core sites, a relationship between stable-isotopes and accumulation rate is assumed based on the saturation vapor pressure (e.g., Petit et al., 1999).This relationship is not well suited to West and coastal Antarctica.For example, in the modern climate, WAIS Divide is 2 ‰ more depleted in δ 18 O than Byrd despite a ∼ 50 % higher accumulation rate.At both Taylor Dome and Law Dome, the Holocene accumulation rates show little resemblance to the stable-isotope records (Monnin et al., 2004;van Ommen et al., 2004).These observations suggest that both spatial and temporal relationships between stable isotopes and accumulation rate are complex.
Ice-flow models forced by accumulation rates inferred from stable isotopes do not exactly match tie points (Ruth et al., 2007;Parrenin et al., 2007) and are commonly adjusted to better match the depth of age markers (e.g., Dreyfus et al., 2007).The Bayesian statistical approach, known as DATICE, (Lemieux-Dudon et al., 2010;Veres et al., 2013;Bazin et al., 2013) allows the thinning function, accumulation rate, and lock-in depth to vary within a tolerance to better reconcile the age of tie points with the initial modeled timescale.We compare ALT with three different published timescales for EDML and discuss the consistency of the interpolations in a later section.
ACCUM and ALT have not previously been used to derive timescales for ice cores.We choose to find the smoothest histories because we want the inferred histories to have the minimum structure required to fit the data.ACCUM minimizes the variability in the inferred accumulation history while ALT minimizes layer-thickness variations, which are related to the accumulation rate by the thinning function.Because the thinning function, which gives the cumulative amount of thinning a layer has experienced, is smooth in time, this leads to both ACCUM and ALT constraining the accumulation-rate variability.ACCUM requires an ice-flow model to estimate the thinning function.Although we use a one-dimensional ice-flow model (Dansgaard and Johnsen, 1969;WAIS Divide Project Members, 2013), a coupled thermo-mechanical model could also be used.We assume that all inputs necessary for the ice-flow model except the accumulation rate are known.The influence of the assumed iceflow model inputs are discussed in a following sub-section.Assumptions about model parameters and ice-flow history are not needed when using ALT.
The underlying assumption that accumulation rate varies smoothly breaks down in situations of abrupt changes, as observed in Greenland at Dansgaard-Oeschger events (e.g., Alley et al., 1997).Even relatively small changes such as the increase in accumulation rate evident in WDC06A-7 between 12.0 and 11.6 kyr are difficult to match without closely spaced tie points; this event is discussed in more detail below.

Smoothest accumulation rate and annual-layer thickness interpolations
We use a standard inverse procedure (Aster et al., 2005) to minimize the data misfit and the smoothness of either the annual-layer thickness or accumulation-rate history for a given trade-off parameter ν, and define a performance index: where m is the parameter being solved for (either accumulation-rate history or annual-layer thickness profile), G(m) is a function that relates m to the depth-age scale, d are the depths of tie points, σ is the standard deviation of the measurement error, which is assumed to be normally distributed, and L is the matrix second derivative operator for calculating smoothness.Bold capital letters denote matrices, bold lower-case letters denote vectors, and || || 2 2 indicates an L2 norm (sum of squares, e.g., Aster et al., 2005).The uncertainty, σ , is a scalar because we assume that all tie-points are equally well known.
The forward problem in this solution procedure, denoted as G(m), can be any function that maps either accumulation rate or annual-layer thickness to a depth-age relationship.In the case of ACCUM, m is the set of accumulation rates at the specified calculation times.The forward problem is a Dansgaard-Johnsen (1969) ice-flow model modified following Dahl-Jensen et al. (2003) to account for basal melting and sliding.In addition to the accumulation rate, the ice-flow model requires histories of basal melting, the fraction of surface motion due to basal sliding, ice thickness, and the kink-height, which dictates the shape of the verticalvelocity profile.We choose this model because it has been applied widely to ice cores and is computationally inexpensive, but any ice-flow model could be used.In the case of the smoothest annual-layer thickness, the function G(m) is simply an integration of the annual-layer thickness profile.The ACCUM interpolation can be written as a linear inverse problem (Lundin, 2012), which allows the uniqueness of the solution to be proven.However, the linear method is also restricted to sites where the ice sheet is frozen to the bed.By using a generalized non-linear solution procedure, we can apply this technique to sites like WAIS Divide that have experienced basal melting, with confidence that the solution is unique based on the linear formulation (see Appendix A).Further, the solution technique can be extended to a range of problems that are not linear, such as finding the ice-thickness history that best matches the depth-age tie points given a known accumulation-rate history (e.g., Price et al., 2007).
We find the best-fit models using a Gauss-Newton iterative inverse procedure (Aster et al., 2005, Sect. 10.1;Ganse, 2013).This is a steepest descent solver that uses the linear relationship between the model perturbations and the data residuals (Appendix A).We fit the known depth-age tie points near exactly because this work is focused on the sensitivity of the timescales to the interpolation method.The trade-off parameter, ν, is chosen such that the maximum misfit in depth at any depth-age tie point is less than the minimum annual-layer thickness in the derived timescale -typically between 1 and 10 mm.We do not discuss the application using the smoothness criteria to improve the age estimates of the depth-age tie points; this application is discussed by Lundin (2012).

Comparison of ACCUM with different ice-flow parameters
We assess the importance of the ice-flow model inputs in the ACCUM method by comparing timescales that assume very different inputs.We use WDC06A-7 to define tie points at the approximate ages of methane ties (Table 1) and match these ages near exactly.In the first case (flank), we use reasonable values for WAIS Divide: no thickness change, a Dansgaard-Johnsen kink height of 20 % of the ice thickness appropriate for flank flow, a basal melt rate of 1 cm yr −1 , and sliding over the bed contributing 50 % of the surface velocity.In the second case (divide), we also use no thickness change, but we specify a Dansgaard-Johnsen kink height of 0.7 for divide flow, and no basal melting or sliding (Raymond, 1983;Conway et al., 1999).Inferred accumulation rates are shown in Fig. 2a.The different ice-flow assumptions result in inferred accumulation rates that differ by more than a factor of 4. The divide case has unrealistically high accumulation rates in the glacial period.Figure 2b shows age differences between the two models for ice at given depths.By construction, the age differences are zero at the tie points.The age differences are also relatively small between the tie points, reaching a maximum of 56 years in the middle of the 6200-year gap between the tie points at 17.8 and 24.0 kyr.The small age differences indicate that ACCUM is relatively insensitive to the prescribed ice-flow parameters.Because the thickness of a layer is the product of the accumulation rate at the age of the layer and the thinning function (Appendix B), changes in the thinning function can be compensated by the accumulation-rate history.The insensitivity of ACCUM to the choice of ice-flow parameters is useful when interpolating timescales.

Assessment of interpolation methods using WAIS Divide data
The annually resolved WDC06A-7 provides the first opportunity to directly assess interpolation timescale methods with an ice core from a location with relatively smoothly varying accumulation.The 100-year running average of annuallayer thicknesses is shown in Fig. 3. Significant centennial variability is evident, but apart from the notable exception between 11.6 and 12 kyr, there is no evidence of abrupt climate changes similar to Dansgaard-Oeschger events in the Greenland ice-core records.Between 11.6 and 12.0 kyr the average annual-layer thickness increases by ∼ 40 %, which is interpreted as a change in accumulation rate (WAIS Divide Project Members, 2013).The increase occurs during the peak of Antarctic Isotope Maximum 0 (corresponding to the end of the Younger Dryas) when the stable-isotope values show little change.Whether this accumulation increase is specific to WAIS Divide or extends across much of Antarctica is unclear because of the lack of annually resolved timescales for other ice cores.This event shows that even though the Antarctic climate varies relatively smoothly, events that deviate from our expectations are possible.
The annual-layer thicknesses of the three interpolation methods are also shown in Fig. 3.The linear interpolation causes step changes in the inferred annual-layer thickness at each tie point.At some tie points, the difference in annuallayer thickness on either side of the tie point is large; the annual-layer thickness changes by nearly 100 % from a 3 cm thickness for ages older than 17.8 kyr to a 6 cm thickness for ages younger.The other two interpolation methods yield smoothly varying annual-layer thicknesses that match each other and the measured annual-layer thicknesses reasonably well.
The differences in age between the annually resolved WDC06A-7 and the interpolations are shown in Fig. 4; the yellow shading shows the uncertainty of WDC06A-7 accumulated from the nearest tie point.The uncertainty accumulates quasi-linearly because either too many or too few years may have been systematically identified (see WAIS Divide Project Members, 2013 for specifics of the WDC06A-7 timescale and associated uncertainties).The ALT and AC-CUM interpolations tend to have similar age differences from WDC06A-7.The linear interpolation, however, tends to be older than WDC06A-7.Because constant layer thickness between tie points is assumed with linear interpolation, layers are too thin on the younger side of the interval and too many years are present.After the approximate midpoint, the layers are then too thick and too few years are present.In contrast, in the interval between 24 and 28 kyr, the measured annuallayer thicknesses increase with age and the linear interpolation places too few layers at younger ages.
All interpolation methods do a poor job of matching WDC06A-7 between 8 and 12 kyr.This is not surprising, since interpolation cannot accurately capture abrupt variations between tie points.Interestingly, linear interpolation provides the best match to WDC06A-7 during this interval because large annual-layer thicknesses near 12 kyr reverse the trend of decreasing layer thickness with age.The other interpolation methods underestimate the age because they yield smaller layer thicknesses at the older side of the interval.
In the interval from 15 to 18 kyr, layer thicknesses decrease rapidly with age on WDC06A-7 due to the combination of layer thinning from ice-flow and the glacialinterglacial accumulation-rate change.The linear interpolation has a particularly large age difference because of the large changes in layer thickness during the interval.In the interval 18 to 24 kyr, the annual-layer thickness varies much less.ALT and ACCUM yield layers that are too thick near 18 kyr, with corresponding ages that are too young.
The mismatch between the linear interpolation and WDC06A-7 is ∼ 374 years, with ages that are outside the estimated uncertainty of WDC06A-7 most of the time (Table 2).Both ACCUM and ALT are improvements: age differences between the two are less than 20 years except during the interval 24 to 17.8 kyr.Interestingly, ALT is slightly better than ACCUM in each metric.The maximum mismatch for ALT is 118 years compared to 151 years for ACCUM.ALT also yields fewer ages that are outside the WDC06A-7 uncertainty.The two methods are expected to be similar The slightly better performance of ALT suggests that the greater complexity of ACCUM is not warranted; this may be particularly true for the deepest ice where the thinning function becomes increasingly uncertain and difficult to predict with ice-flow models.ACCUM will likely improve in cases where the ice physics or ice-flow histories are better known; using two or three dimensional ice-flow models constrained by dated internal layers imaged by radar is a promising approach (e.g., Waddington et al., 2007;Steen-Larsen et al., 2010).However, in the following we will use ALT to assess existing ice-core timescales because it is simpler to apply.

Age uncertainty due to interpolation
Comparison of the interpolation methods provides a framework for estimating the interpolation uncertainty for other Antarctic ice-core timescales.In general, the age difference between the interpolation and WDC06A-7 increases farther away from a tie point (Fig. 3).We use the age differences to estimate the rate at which different interpolation methods accumulate age uncertainty and calculate the rate of accumulating uncertainty as the absolute value of the age difference divided by years from closest tie point.For instance, if the age difference is 50 years at a point 500 years from the closest tie point, then the uncertainty has accumulated at a rate of 10 years per hundred years.This is done for all of the ages between 2 and 29.6 kyr.We exclude the topmost interval because the interpolated firn-density profile does not exactly match the actual density profile.We then find the cumulative fraction of ages for each rate of accumulating uncertainty (Fig. 5).The cumulative fractions can be used as rough estimates of the 1-σ and 2-σ rates of accumulating uncertainty.We emphasize that these are not formal statistical uncertainties because they are based on a single timescale; however, they provide a rough estimate in the absence of other quantitative information about timescale uncertainties away from tie points.
For ALT and ACCUM, the rate of accumulating age uncertainty is less than 4 years per hundred 68 % of the time, and less than 10 years per hundred 95 % of the time (Fig. 5).The maximum rate of accumulating uncertainty is about 20 years per hundred and occurs at the abrupt accumulation increase at ∼ 12 kyr.When the annual-layer thickness varies more rapidly, there tend to be higher rates of accumulating uncertainty.Linear interpolation accumulates uncertainty at a significantly greater rate than ACCUM and ALT.The rates of accumulating uncertainty for the linear interpolation are less than 8 years per hundred years 68 % of the time and less than 31 years per hundred years 95 % of the time.Linear interpolation accumulates uncertainty at a rate greater than 50 years per hundred years approximately 3 % of the time.

Duration uncertainty due to interpolation
For many analyses, it is important to know both the age and duration of a climate event.The duration depends on the number of annual layers between the measured depths of the onset and termination of the event.Interpolation tends to underestimate variability in duration because it cannot capture high-frequency variations in annual-layer thickness.Interpolation can also artificially shorten or lengthen the duration of climate events, making them appear more or less abrupt.The effect can be large for linear interpolation when the events occur at the beginning or end of an interval during which the actual annual-layer thicknesses varied (e.g., 15 to 18 kyr at WAIS Divide, Fig. 3).An additional issue with linear interpolation is when a climate event begins before a tie point and ends after.For instance, the annual-layer thicknesses increase from 3 cm for ages older than the 17.8 kyr tie point to 6 cm for ages younger (Fig. 3).For equally spaced measurements in depth, twice as much time is packed between measurements on the older side of the tie point than between measurements on the younger side of the tie point.Therefore, features just older than 17.8 kyr will have twice the duration as features just younger.This may influence detection of times of significant change by artificially altering the rate of change of a climate proxy.

Age bias due to interpolation
The choice of interpolation method can also bias the average age of the timescale.Understanding the magnitude of potential bias is relevant for analyses of the phasing of climate proxies from different ice cores, such as between a composite Antarctic temperature record from many cores and carbon dioxide records measured from individual cores (e.g., Pedro et al., 2012).ACCUM and ALT produce ages that are older than WDC06A-7 about as often as they are younger (Fig. 4); however, the linear interpolation is nearly always older (84 % of the time).This occurs because the layer thickness predominantly decreases with age due to thinning from ice flow.The constant annual-layer thicknesses between tie points resulting from linear interpolation cause the annual layers to be too thin near the younger tie point and therefore put too much time into the younger half of the depth interval.The situation reverses in the older half of the interval, but these layers are still too old at each depth, recovering the correct age only at the bottom of the interval.Therefore, the linear interpolation is too old for the entire interval.For WDC06A-7, the average bias is 55 years for the timescale as a whole (Table 2).
The bias due to interpolation tends to be larger for time periods that contain fewer tie points because there are fewer opportunities for older and younger age differences to cancel out (Table 2).The bias is greatest for linear interpolation during periods with large variations in annual-layer thicknesses.For instance, the increasing annual-layer thicknesses of the deglacial transition (21 to 12 kyr) result in ages biased 181 years too old for linear interpolation compared to WDC06A-7.In contrast, the ACCUM bias increases to three decades and the ALT bias is negligible (Table 2).

Comparison of ALT with timescales for EDML
We compare three published timescales for EDML with an interpolation using ALT.The three timescales are AICC2012 (Veres et al., 2013); LD2010 (Lemieux-Dudon et al., 2010);  (Veres et al., 2013), LD2010 (Lemieux-Dudon et al., 2010), and EDML1 (Ruth et al., 2007) -and the ALT interpolation based on ice ages from AICC2012.(B) Age differences between AICC2012 and ALT interpolation.Age differences between other timescales are not shown because only the AICC2012 and ALT use the same ages for tie points.Dashed vertical lines are tie-point ages.and EDML1 (Ruth et al., 2007).To derive the ALT timescale, we use the ice ages at times of distinctive methane features from the AICC2012 timescale (Table 3) as well as a few tie points in the Holocene.The annual-layer thickness profile from AICC2012 and ALT are shown in Fig. 6a.For most of the timescale, annual-layer thicknesses agree closely, and the age differences are less than 100 years (Fig. 6b).The largest difference occurs between 12.8 and 9.2 kyr, when ALT does not produce thick annual layers around 11.5 kyr.The thick layers in the AICC2012 are driven by high isotope values at the Antarctic Isotope Maximum 0 peak.
Annual-layer thickness profiles from the EDML1 and LD2010 timescales are also shown in Fig. 6a.These two timescales use the same abrupt methane variations as AICC2012 for ages older than 11.5 kyr (ice age), but ages younger in AICC2012 are based on volcanic sulfur matches with the North Greenland ice core and annual-layer counting (Veres et al., 2013).The inferred layer thicknesses using ALT are, on average, more similar to those of AICC2012 than are the annual-layer thicknesses of either EDML1 or LD2010 to AICC2012.EDML1 predicts thicker layers at both Antarctic Isotope Maximums 0 (∼ 12 ka) and 1 (∼ 15 ka).LD2010 predicts a 25 % increase in annual-layer thickness between 8 and 7.2 kyr when the other timescales show little change.The smoothness requirement causes ALT to miss high-frequency structure in layer-thickness profiles, but that requirement also avoids creating large variations in annual-layer thickness that may not be real.That the differences in layer thickness among the three EDML timescales are larger than the differences between ALT and AICC2012, indicating ALT is a reasonable interpolation for EDML.Age differences between ALT and EDML1 and LD2010 are not shown because the age of the tie points do not match those of AICC2012; thus an age comparison is not informative about the interpolation method.The similarity of the annual-layer thickness profiles among the four timescales suggests that both the Bayesian inverse and ALT methods produce robust timescales for EDML.

Application to the Byrd, Siple Dome, and Law Dome ice-core timescales
Comparison of interpolation methods with the annually resolved WDC06A-7 timescale indicates that linear interpolation can lead to large age differences.ALT and ACCUM performed substantially better.Here we assess three timescales that are based on near-linear interpolations: Byrd, Siple Dome, and Law Dome.We use ALT because it does not require any assumptions about ice flow.Ages of tie points for each timescale are derived from the original timescales either directly from the stated tie points or indirectly from the abrupt variations in annual-layer thickness (Table 4).We do not consider Taylor Dome because of the large uncertainty of the tie points (Mulvaney et al., 2000).For Taylor Dome, improving only the interpolation is unlikely to yield insight into the climate history.
Annual-layer thickness profiles are shown in Fig. 7a, c, and e for Byrd, Siple Dome, and Law Dome, respectively.All of the timescales have abrupt changes in annual-layer thickness at tie points, similar to the linear interpolation of WAIS Divide presented in Fig. 3, but also have structure due to the specifics of the timescale construction.The Byrd gas timescale was first developed based on a Monte Carlo method that maximized the correlation between methane records (Blunier and Brook, 2001).The ice timescale was then derived by adding the age, calculated from accumulation and temperature estimates based on the stable-isotope values.The abrupt variations in methane drive the correlation, such that the resulting timescales have similarities with linearly interpolated timescales.The high-frequency variations in annual-layer thicknesses in the ice timescale result from both the correlation of methane records and the time-varying age values.The Siple Dome gas timescale was derived by matching the methane record to the Greenland methane composite record at distinct tie points then linearly interpolating between the tie points (Brook et al., 2005).High-frequency variations in annual-layer thickness result from different age values.The Law Dome timescale was developed using a Dansgaard-Johnsen (1969) ice-flow model to calculate layer thinning and assumed constant accumulation rates between the tie points (van Ommen et al., 2004;Pedro et al., 2011).Because there is little variation in the layer thinning in the deep ice, the resulting timescale is similar to linear interpolation (see Fig. 7e).
We have re-interpolated the three timescales with ALT.While most of the tie points in the original timescales were used, Byrd and Siple Dome required a few exceptions: we used only one, instead of two, tie points at 24, 27.5 and 29 kyr (Dansgaard-Oeschger events 2, 3 and 4) to avoid inferring large, unrealistic variations in layer thickness over short intervals.
ALT annual-layer thickness profiles are shown in red in Fig. 7a, c, and e. Age differences between the original nearlinear interpolations and ALT are shown in black.All three timescales show significant age differences.The largest differences are in the Siple Dome timescale; ages shift by as much as 1200 years around 20 kyr, which is within the estimated 2000 years uncertainty for the timing of the abrupt isotope change ∼ 22 kyr (Brook et al., 2005).For Byrd, age differences are up to 1000 years at about 24 kyr.At both Siple Dome and Byrd, changes in the timescales more recent than The stable-isotope records on the ALT timescales are shifted down by 1 ‰.
18 kyr are less than 200 years.At Law Dome, the timescale ends at 21 kyr.There is no tie point at the onset of the deglacial rise (∼ 18 kyr), which results in a relatively long span, 16.2 to 19.5 kyr, when annual-layer thickness decreases by a factor of 2 between the tie points.The result is a large age discrepancy of over 500 years centered at 18.1 kyr.As shown with the WDC06A-7 timescale, linear or near-linear interpolation performs poorly when the annual-layer thicknesses change significantly between tie points.Stable-isotope records are shown in Fig. 7b, d, and f on the original (blue) and ALT (red) timescales.For Byrd, the most significant change is that Antarctic Isotope Maximum 2 (∼ 24 kyr) becomes more compact and similar in character of the WAIS Divide and EDML records (EPICA Members, 2006;WAIS Divide Project Members, 2013).Changes in the glacial-interglacial transition part of the record are smaller, although ages during the Antarctic Cold Reversal are shifted up to 150 years younger.For Siple Dome, the 3 ‰ δ 18 O increase is shifted younger by 1200 years from 21.8 to 20.6 kyr.As with Byrd, the changes during the glacial-interglacial transition are smaller but up to 250 years younger at 17.2 kyr.For Law Dome, the largest shift is centered at 18.1 kyr.The onset of deglacial warming at Law Dome is examined in more detail below.
Age differences between the near-linear interpolation and ALT indicate that linear interpolation consistently yields older ages during the deglacial transition.Except for a few short periods, the ALT timescales are younger between 21 and 12 kyr.The bias towards older ages using a linear interpolation occurs because annual layers are too thin on the younger side of an interval during times of decreasing (in age) accumulation rates.On average, the linear timescales are too old by 40 years for Byrd, 240 years for Siple Dome, and 150 years for Law Dome.
The magnitude of the linear interpolation biases are similar to the magnitude of the century-scale lead (−56 to 381 years) of Antarctic temperature to carbon dioxide found by Pedro et al. (2012) and suggest interpolation biases have the potential to affect the details of such analyses.Whether there is any significant impact on the Pedro et al. (2012) conclusion is unclear because timescales for both the ice (Antarctic temperature) and gas (carbon dioxide) were constructed with near-linear interpolation: Byrd, Siple Dome, and Law Dome are three of the five δ 18 O records in the Antarctic temperature stack and the two gas timescales for carbon dioxide are Byrd and Siple Dome.Redoing the lag analysis is beyond the scope of this paper, but the potential biases indicate future investigations of century-scale leads and lags should include interpolation uncertainties.
We cannot verify whether the ages produced by ALT are more accurate than the original timescales but the ALT timescales are valid alternatives.The large differences in the timescales due to the chosen interpolation method highlight the increased uncertainty in the timescale away from the tie points.Caution is needed when interpreting climate changes based on interpolated timescales.

Example: onset of the glacial-interglacial transition at Law Dome
The onset of deglacial warming at Law Dome is marked by a rapid increase of 3 ‰ δ 18 O in ∼ 500 years.The timing of the onset was determined to be 17.84 ± 0.32 kyr (Pedro et al., 2011) on the original timescale (blue in Fig. 8); the blue vertical line marks the onset and the blue shading shows the stated uncertainty.Figure 8 also shows the same stableisotope record on the ALT timescale.The onset of deglacial warming is shifted 550 years younger to 17.29 kyr.The large age difference occurs because this climate feature is nearly midway between tie points, and the annual-layer thickness changes by a factor of 2 throughout the interval.At this climate feature, the interpolation uncertainty exceeds the stated uncertainty of the age and correlation uncertainties at the tie points.Both the timing of the climate feature and stated uncertainty should be revised.Assuming that the rate of accumulating uncertainty for linear interpolation developed with the WAIS Divide timescale (31 years per hundred, Fig. 5) is appropriate for the Law Dome timescale, the total uncertainty for this climate feature is 830 years: the 320 years stated uncertainty plus an interpolation uncertainty of 510 (0.31 × 1640) years from the closest tie point.2011) timescale (blue) and on ALT timescale (red and shifted down by 1 ‰).Black vertical dashed line is the depth-age tie point.Blue vertical line and shading are the respective timing and uncertainty of the onset of deglacial warming defined by Pedro et al. (2011).The red vertical line marks the timing of the onset of deglacial warming with ALT.The blue horizontal bar is increased uncertainty by adding an interpolation uncertainty appropriate for linear interpolation.The red horizontal bar is the total uncertainty for the ALT timescale.See text for description of uncertainty calculations.The two circled areas show the different durations of climate events that can result from using different interpolation methods.The duration of the cooling with ALT is about half its duration with linear interpolation.
In the ALT timescale, the onset of deglacial warming begins at 17.29 kyr and the total uncertainty would be 430 years, the 320-year stated uncertainty and 110-year (0.1 × 1100) interpolation uncertainty.The distance to the closest tie point changes with the interpolation method, which affects the magnitude of the interpolation uncertainty.This reinforces that the uncertainties are estimates and should not be interpreted as a precise quantification of the total age uncertainty.
A second effect is that the duration of events is different among the interpolation schemes.This is illustrated by the circled data in Fig. 8.The duration of the cooling beginning at 16.6 kyr on the linear timescale is nearly twice as long as on the ALT timescale due to the large difference in annuallayer thickness between the two timescales (Fig. 7e).Statistical tests identifying times of significant change (e.g., Pedro et al., 2011;WAIS Divide Project Members, 2013) may be affected when the annual-layer thicknesses change quickly at tie points.

Conclusions
Three methods for interpolating between age markers were tested using the annually resolved timescale for the WAIS Divide ice core (WDC06A-7).Linear interpolation results in unrealistic variations in annual-layer thicknesses at tie points, causing both large age errors and abrupt apparent variations in the duration of climate events.In addition, the linear interpolation resulted in a bias to older ages.The mismatch improved using interpolations based on either the smoothest annual-layer thickness (ALT) or smoothest accumulationrate history (ACCUM).ALT performed slightly better than ACCUM but both give sufficiently similar results, such that the choice of which to use will depend on the application.If only a timescale is desired, ALT is simpler to implement, but if accumulation-rate estimates are desired, for instance to help constrain the age, then ACCUM is required.The rate of accumulating uncertainty between tie points using ALT was 10 years per hundred years, compared to 31 years per hundred years for linear interpolation.
The existing Byrd, Siple Dome, and Law Dome timescales were based on near-linear interpolation.Re-interpolation of these timescales using ALT and similar tie points produced timescales that were on average younger during the glacialinterglacial transition.For Siple Dome, the abrupt 3 ‰ increase in δ 18 O at ∼ 22 kyr (Taylor et al., 2004) was shifted 1200 years younger; while the timing of this event has significant uncertainty due to the accuracy of the methane tie points, the large shift in age due to interpolation indicates that caution should be exercised when using the Siple Dome record in analyses of the onset of Antarctic deglacial warming.At ages between tie points, the interpolation uncertainty can exceed uncertainty at the tie points.For Law Dome, the timing of the onset of deglacial warming is shifted by 540 years, indicating that the interpolation uncertainty is nearly double the stated 320 years which was determined from the tie point uncertainty (Pedro et al., 2011).Uncertainty in interpolation is often omitted in analyses; it should be considered when determining the timing and duration of climate features in Antarctic ice-core records.
The Supplement related to this article is available online at doi:10.5194/cp-10-1195-2014-supplement.
ACCUM minimizes the second derivative of the accumulation rate: ALT minimizes the second derivative of the annual-layer thickness: Equation ( B3) can be rewritten after substituting in Eq. ( B1) Unfortunately, ALT does not compute a thinning function or an accumulation history needed to evaluate the three terms on the right-hand side of Eq. (B4).However, ACCUM does calculate all of these terms and we can use these to better understand how the ALT minimization differs from the AC-CUM minimization.Figure B1 shows relative values for the three terms on the right-hand side of Eq. (B4), evaluated for ACCUM as a function of age on the WDC06A-7 timescale.The third term, which includes the second derivative of the accumulation rate, accounts for 70 % of the curvature of λ(A) and explains why ALT and ACCUM yield similar, but not identical, interpolations.The second derivative of the accumulation rate is multiplied by the thinning function; because the thinning function is unity at the surface, ALT emphasizes minimizing the curvature of the younger accumulation rates relative to ACCUM.
The first and second terms each account for 15 % of the total.The first term is the second derivative of the thinning function weighted by the accumulation rate.Hence ALT will minimize the thinning function with added emphasis when the accumulation rate is higher, such as in the mid-Holocene for WAIS Divide.The second derivative of the thinning function is also indirectly controlled by variations in accumulation rate because a major cause of variations in the thinning function is the change in accumulation rate.The second term, with the two first derivatives, has a similar structure to the first term.The largest age difference between ALT and AC-CUM is in the 18 to 24 kyr interval.Interestingly, this interval has near-zero values for both the first and second terms, making it difficult to explain why the ALT and ACCUM differ during this period.The difference is likely related to weighting of the second derivative of the accumulation rate by the thinning function in the third term, which lessens the smoothness constraint during this period in ALT.
It is not clear that ALT has a predictive advantage to AC-CUM, or whether the slightly better performance at WAIS Divide is due to chance.There is no obvious reason why minimizing the curvature of the thinning function would lead to improved predictive capability.While this is an area that could be further explored, the small differences between the ALT and ACCUM methods indicate that there will be little improvement in the final timescale.The choice of whether to use ALT or ACCUM will most likely depend on the application.If only a timescale is desired, ALT is simpler to implement.If accumulation rate estimates are desired, for instance to help constrain the age, then ACCUM is required.

Figure 2 .
Figure 2. Accumulation rate (top panel) inferred with the ACCUM method with ice flow parameters for either a flank or divide site.The age difference (bottom panel) between the two interpolations.

Figure 3 .
Figure 3.A 100-year running average of annual-layer thicknesses for the WDC06A-7 compared with inferred annual-layer thickness for the different interpolations.Black vertical dashed lines are tie points.Inset in upper left shows detailed view for 29.6 and 17 kyr intervals.

Figure 4 .
Figure 4. Age differences between the interpolation methods and the WDC06A-7 timescale.Yellow shading is the accumulated uncertainty of WDC06A-7 from the closest tie point (WAIS Divide Project Members, 2013).
Figure 6.(A)Annual-layer thicknesses for three EDML timescales -AICC2012(Veres et al., 2013), LD2010(Lemieux-Dudon et al., 2010), and EDML1(Ruth et al., 2007) -and the ALT interpolation based on ice ages from AICC2012.(B) Age differences between AICC2012 and ALT interpolation.Age differences between other timescales are not shown because only the AICC2012 and ALT use the same ages for tie points.Dashed vertical lines are tie-point ages.

Figure 7 .
Figure 7. Annual-layer thicknesses of original timescales (blue) and ALT interpolation (red) in (A), (C), and (E).Age differences are shown in black on the original timescale.The stable-isotope records on the original (blue) and ALT (red) timescales are shown in (B), (D), and (F).The stable-isotope records on the ALT timescales are shifted down by 1 ‰.

Figure 8 .
Figure 8. Detail of Law Dome δ 18 O record on the Pedro et al. (2011) timescale (blue) and on ALT timescale (red and shifted down by 1 ‰).Black vertical dashed line is the depth-age tie point.Blue vertical line and shading are the respective timing and uncertainty of the onset of deglacial warming defined byPedro et al. (2011).The red vertical line marks the timing of the onset of deglacial warming with ALT.The blue horizontal bar is increased uncertainty by adding an interpolation uncertainty appropriate for linear interpolation.The red horizontal bar is the total uncertainty for the ALT timescale.See text for description of uncertainty calculations.The two circled areas show the different durations of climate events that can result from using different interpolation methods.The duration of the cooling with ALT is about half its duration with linear interpolation.

Figure B1 .
Figure B1.The three terms from the right side of Eq. (B4) evaluated for the ACCUM interpolation of WAIS Divide WDC06A-7.

Table 2 .
Interpolation assessment for WDC06A-7.This is the percentage of interpolation ages that differ from the WDC06A-7 ages by more than the WDC06A-7 uncertainty at that age.since the ACCUM interpolation produces a smooth annuallayer thickness.Reasons why ALT and ACCUM might differ are discussed in Appendix B.
We define the rate of accumulating uncertainty as the absolute value of the age difference between an interpolation and the WDC06A-7 timescale at each WDC06A-7 age divided by the age difference from the nearest tie point.The cumulative fraction is the fraction of the total data set for which the rate of accumulat- ing uncertainty is below the threshold value on the x axis; this is similar to a cumulative distribution function.The horizontal lines indicate 0.68 (1σ ) and 0.95 (2σ ) levels.