A first chronology for the North Greenland Eemian Ice Drilling ( NEEM ) ice core

A stratigraphy-based chronology for the North Greenland Eemian Ice Drilling (NEEM) ice core has been derived by transferring the annual layer counted Greenland Ice Core Chronology 2005 (GICC05) and its model extension (GICC05modelext) from the NGRIP core to the NEEM core using 787 match points of mainly volcanic origin identified in the electrical conductivity measurement (ECM) and dielectrical profiling (DEP) records. Tephra horizons found in both the NEEM and NGRIP ice cores are used to test the matching based on ECM and DEP and provide five additional horizons used for the timescale transfer. A thinning function reflecting the accumulated strain along the core has been determined using a Dansgaard– Johnsen flow model and an isotope-dependent accumulation rate parameterization. Flow parameters are determined from Monte Carlo analysis constrained by the observed depth-age horizons. In order to construct a chronology for the gas phase, the ice age–gas age difference ( 1age) has been reconstructed using a coupled firn densification-heat diffusion model. Temperature and accumulation inputs to the 1age model, initially derived from the water isotope proxies, have been adjusted to optimize the fit to timing constraints fromδ15N of nitrogen and high-resolution methane data during the abrupt onset of Greenland interstadials. The ice and gas chronologies and the corresponding thinning function represent the first chronology for the NEEM core, named GICC05modelext-NEEM-1. Based on both the flow and firn modelling results, the accumulation history for the NEEM site has been reconstructed. Together, the timescale and accumulation reconstruction provide the necessary basis for further analysis of the records from NEEM.


Introduction
The 2540 m long North Greenland Eemian Ice Drilling (NEEM) ice core was drilled during 2008-2010 through the ice sheet in Northwest Greenland (77.45 • N, 51.07 • W, surface elevation 2479 m).As for all palaeoclimate sediment records, a (depth, age) relationship, or timescale, is needed to fully utilize and interpret the measurements performed on NEEM ice and gas samples.

S. O. Rasmussen et al.: A first chronology for the NEEM ice core
The present-day annual snow accumulation at the NEEM site is 0.22 m ice equivalent and, therefore, the core has annual layering similar to the GRIP and GISP2 ice cores.Once high-resolution measurements from the NEEM core become available, annual layer identification in the NEEM core may be undertaken either by the manual identification of annual layers as performed on data from the DYE-3, GRIP, and NGRIP cores in the construction of the Greenland Ice Core Chronology 2005, hereafter GICC05, (Vinther et al., 2006;Rasmussen et al., 2006), or by automated methods (Winstrup et al., 2012).
At this point, however, data for annual layer identification are not available.Thus, to provide a first timescale for the NEEM core, the GICC05 timescale and its model extension GICC05modelext is adopted as is from the NGRIP core using linear interpolation between common reference horizons, or match points, as identified from the electrical conductivity measurement (ECM), dielectrical profiling (DEP), and tephra records of the two cores.
Analysis of several ice-and gas-based records from the NEEM core shows that the lower parts of the core, which contain ice from the Eemian interglacial, have been subject to folding, but no such discontinuities have been identified in the upper 2203.6 m of the core (NEEM community members, 2013).The GICC05modelext timescale only extends to about 122 ka b2k (thousand years before 2000 CE, Wolff et al., 2010) -corresponding to the lower limit of the NGRIP core -and therefore does not cover the entire NEEM Eemian record.The lower section of the NEEM core below 2203.6 m, having complex stratigraphy, was therefore tied to the Antarctic EDML1 timescale by NEEM community members (2013) using the methane (CH 4 ) and δ 18 O air records (Capron et al., 2010).
The work presented here focuses entirely on the stratigraphically uncompromised section of the NEEM core, where it is reasonable to assume the existence of an unambiguous and continuous mapping between NEEM and NGRIP depths, which allows us to align common patterns of peaks in the ECM and DEP records of the two cores.

Data
This work is based on the transfer of the existing GICC05modelext timescale, and on previously unpublished ECM and DEP data sets as well as a number of tephra horizons.The data used for the age calculations in Sect.3.4 are presented by Chappellaz et al. (2013) and Guillevic et al. (2013).

The GICC05 and GICC05modelext timescales
Annual layers can be identified in all deep ice cores from Greenland, and due to differences in accumulation rates, ice thickness, temperature, ice flow conditions, and availabil-ity of data in sufficient resolution, different deep ice core records are best suited for annual layer identification in different age or depth ranges.For example, with a present-day annual accumulation of 0.56 m of ice equivalent, the DYE-3 ice core is optimal for dating the most recent millennia as the DYE-3 δ 18 O record contains a well-developed annual signal reflecting the seasonal temperature variation.However, the high accumulation rates lead to faster thinning of layers with depth (Vinther et al., 2006).In contrast to DYE-3, the present-day annual accumulation at the NGRIP site of 0.19 m of ice equivalent means that the NGRIP δ 18 O record only marginally resolves the annual layers even at the top of the core.However, the combination of the availability of highresolution continuous flow analysis (CFA) impurity records and comparably low thinning rates in the glacial section due to bottom melting mean that the annual layer thickness remains above 1 cm allowing annual layer identification back to approximately 60 ka b2k (Svensson et al., 2008).
The idea behind the design of the GICC05 timescale is to exploit these differences, creating a multi-core timescale based on layer counting, using those cores and records that optimally resolve the annual layers at a given time period.The resulting timescale is based on identification of annual layers in δ 18 O records in the top section, using parallel data from DYE-3 (most recent 8.3 ka), GRIP (most recent 3.8 ka), and NGRIP (most recent 1.8 ka) linked together using common volcanic signatures (Vinther et al., 2006).Between 11.7 ka b2k (the onset of Holocene) and 7.9 ka b2k, GRIP high-resolution CFA records of ammonium, hydrogen peroxide, and calcium ions were used for annual layer counting supplemented by short sections of GRIP δ 18 O data.In the section older than 10.3 ka b2k, the multi-component NGRIP data set is also available (Rasmussen et al., 2006).In the Holocene many different data series have been used and overlaps between the different data series ensure maximum consistency.
Glacial annual layer thicknesses are much smaller than Holocene values due to low glacial accumulation rates and flow-induced thinning and therefore only the NGRIP CFA, ECM and visual stratigraphy (VS) data can be used for annual layer identification in the glacial (Andersen et al., 2006b;Svensson et al., 2006Svensson et al., , 2008)).Beyond 60 ka b2k, the counted GICC05 timescale is extended with the modelled ss09sea06bm timescale (North Greenland Ice Core Project members, 2004) shifted 705 yr towards younger ages to ensure continuity (Wolff et al., 2010).The combined timescale has later been named GICC05modelext, and annual layer counting across Greenland Stadial 22 (GS-22) using new high-resolution CFA data largely confirms the validity of the GICC05modelext estimates of the duration of GS-22 (Vallelonga et al., 2012).
With the aim of providing a common chronological framework for the Greenland cores, the GICC05 timescale has been transferred to the Renland and Agassiz ice cores (Vinther et al., 2008) and to the GRIP and GISP2 ice cores in the interval 32.45-11.7 ka b2k (Rasmussen et al., 2008) and further extended to 48 ka b2k for GRIP in Blockley et al. (2012) using an approach similar to the one applied here.Ongoing work will extend the GICC05modelext to the full length of the stratigraphically uncompromised sections of the GRIP and GISP2 cores, and GICC05 ages are also used in integrated inverse modelling efforts aiming at constructing a common chronological framework for Antarctic and Greenland ice cores (Lemieux-Dudon et al., 2010;Veres et al., 2012) Once more data become available, NEEM data can potentially contribute to the refinement of the GICC05modelext timescale, especially in the section older than approximately 7.5 ka b2k.Here the resolution of the DYE-3 δ 18 O record becomes increasingly marginal for resolving annual layers, and even when the GRIP CFA data becomes available in the section older than 7.9 ka, the robustness of the annual layer count is adversely influenced by the fact that the GRIP CFA set consists of only two independent parameters at moderate resolution.The section from about 10.3 ka b2k (deeper than which NGRIP CFA measurements are available) to 7.5 ka b2k is the weakest section data-wise, and although this time period is partially overlapping with the brittle zone in the NEEM core, there is potential for a strengthening of the annual layer identification across this interval once NEEM data are available.Another potential addition to the GICC framework based on NEEM data is semi-automatic and objective annual layer identification using VS and CFA data, possibly in parallel between the NEEM and NGRIP cores, for example, using the method developed by Winstrup et al. (2012).
The uncertainty in GICC05 annual layer identification is given as the so-called maximum counting error (MCE) derived from the number of uncertain annual layers, each counted as 1/2 ± 1/2 yr.Comparisons with other palaeoclimate records indicate that the GICC05 has little or no systematic bias and confirms that the MCE probably is a conservative estimate of the total age uncertainty of the timescale (Svensson et al., 2008).In this work, we follow the same notation and report all GICC05 ages with corresponding MCE values.It should be stressed that the MCE is not a Gaussian uncertainty measure and that MCE values strictly speaking cannot be compared to Gaussian-style dating uncertainties of other records.However, in the absence of a more appropriate uncertainty estimate, we recommend that the MCE is regarded as the 2σ uncertainty of GICC05 in cases where Gaussian uncertainties are needed (Andersen et al., 2006b).

Electrical conductivity measurements (ECM)
Electrical conductivity measurements (ECM) were conducted in the field during the NGRIP (Dahl-Jensen et al., 2002) and NEEM ice core campaigns.For both cores, the ECM set-up and procedure are similar to the one described in Hammer (1980): first, a flat surface along the depth axis of the ice core is polished with a microtome knife to produce a level and freshly cleaned surface.The hand-held ECM instrument consists of a pair of brass electrodes spaced ∼ 1 cm apart, between which a DC voltage difference of 1250 V is maintained.This device is moved along the ice in a steady fashion, while maintaining a constant pressure of the electrodes onto the ice surface.The electrical field causes H + ions in the ice to be displaced resulting in a transient displacement current which is measured.Steady movement of the ECM instrument over the ice surface prevents the ice from being polarized leading to a reduced current over time (Hammer, 1980).Based on comparisons to high-resolution chemistry measurements, it has been established that the recorded ECM signal is almost entirely a response to the acidic components of the ice, even in the presence of large concentrations of neutral salt (Moore et al., 1992).Hence, ECM is a fast and non-destructive technique for providing a high-resolution acidity profile for the ice core.
The ECM technique is often used as one of the primary tools for establishing an event stratigraphy for volcanic eruptions, which are detectable as high-acidity layers in ice cores.ECM was recorded continuously along the NGRIP and NEEM ice cores.For NGRIP, the upper part (down to 1372 m) of the ECM profile was measured on the NGRIP1 core, whereas the lower part (1346-3085 m) is based on data from the NGRIP2 core.Inaccuracies in the logging of the NGRIP cores have given rise to a small offset (∼ 0.43 m) between corresponding volcanic events recorded in the two cores in the overlapping section, with NGRIP1 data being assigned the greater depth (Hvidberg et al., 2002).

Calibration
The ECM signal is temperature dependent and was calibrated to a standard temperature of −14 • C by using an Arrhenius law with an activation energy of 0.23 eV (Neftel et al., 1985).
The exact conversion between ECM and ice acidity is ambiguous.Even for a specific set of electrodes, the data display a large amount of scatter, and several calibration curves have been suggested based on measured H + ion concentrations.The initial calibration curve suggested by Hammer (1980) uses a square-law relationship, but the appropriate parameters in this calibration seem to be ice-core dependent (Moore et al., 1992).This is partly a result of different ambient conditions and operational differences (e.g.pressure applied to the electrodes) during measurements.Also, the exact procedures concerning storage and cleaning of the ice core surface may play a role.For the NGRIP and NEEM cores, many operators have been in charge of the ECM measurements, and the ice core was subjected to varying conditions.Hence, the calibration is likely to be variable with depth.Here, the ECM current i (in µA) has been converted to ice acidity (in µequiv.H + kg −1 ) using the relationship [H + ] = 0.045 × i 1.73 as suggested by Hammer (1980).For our purposes, however, a further fine-tuning of the ECM calibration is unnecessary: the synchronization of the two cores does not rely on absolute values, but solely on the recognition of similar patterns of peaks in the two ice core acidity records.

Data quality
The quality of the ECM measurements is generally higher in ice than in the firn, where the lower density makes it difficult to maintain constant contact area between ice crystals and electrodes.This leads to generally lower currents and enhanced scatter in data from the uppermost part of the ice core.
To ensure the quality of the ECM data, the ECM profile was measured twice for each core section and the two profiles were subsequently checked for agreement.In case of major dissimilarities, measurements were repeated until reproducible results were obtained and a representative run could be selected.In this way, the reproducibility of the ECM data was continuously checked.
The manually operated ECM device is fast and simple, and it allows the scientific investigator to take advantage of the best-preserved section of the ice core surface by manoeuvring the electrodes around cracks or other areas of poorquality ice.Breaks across the entire core diameter cannot be avoided, and can be recognized as abrupt dips in the ECM profile.The location of breaks was logged as part of the ECM measurements, and data from these sections have subsequently been removed from the resulting ECM data sets.
A disadvantage of the handheld electrode approach is that the depth assignment is only accurate to 1-2 cm, which is a consequence of the flexibility in the electrode set-up together with the uncertainty of the registration of the move-ment of the electrodes along the length of the core.The depth uncertainty of 1-2 cm of the ECM data makes a detailed comparison with, for example, high-resolution chemistry measurements of the ice core challenging.Data have subsequently been interpolated to a depth resolution of 1 cm for NGRIP1 and NEEM data and 1 mm for NGRIP2 data, with the different resolution reflecting changes in the data post-measurement treatment.
Despite the rather low accuracy of the depth assignment, the high resolution of the ECM records is important for the synchronization of the ice cores, as it permits robust pattern recognition of features found in both cores.

Dielectrical profiling (DEP)
The dielectrical stratigraphy of the NEEM and NGRIP cores were recorded in the field with the dielectric profiling (DEP) technique (Wilhelms et al., 1998).Comparisons with highresolution chemistry measurements suggest that the DEP signal in Greenland mainly reflects the acidity of the ice and the core's content of ammonium and salts (Moore et al., 1994).Ice core sections of 1.65-2.2m were usually analysed no more than one day after drilling.The brittle zone ice could not be handled directly after drilling and was stored and processed during the subsequent season, but all the data used here comes from below the brittle zone.
At NEEM, approximately every day a free air measurement was recorded for data calibration.For each measurement the cable stray capacitance and conductance -as determined in the corresponding daily free air measurements -are subtracted from the respective capacitance and conductance records.These measurements are not available for NGRIP, and the missing correction for cable stray conductance means that the NGRIP data set still contains a periodic variation with amplitude of up to 2 % of the baseline conductivity and a wavelength of 1.65 m.After multiplying the conductance with the vacuum permittivity, permittivity and conductivity are calculated by division with the free air capacitance.
The relative permittivity increases during the initial densification in the upper part of the core and remains just above three in deeper core sections.Therefore, permittivity values significantly below this threshold are good indicators of sections with bad core quality influenced, for example, by breaks or missing pieces.Data from these sections are not trustworthy, and corrected data sets are created by removing sections where the permittivity drops below a chosen threshold (see Fig. 1).
The NEEM and NGRIP DEP data sets are provided in processed form without correction for the temperature variation in the science trench or in the core-buffer.Therefore, care should be taken when interpreting the absolute levels between sections.For the identification of patterns of peaks between different cores, however, this limitation is not important.

Tephra horizons
Glacial ice from both the NEEM and NGRIP ice cores is currently being intensively investigated with the aim of developing a high-resolution tephrochronological framework, linking the cores together and providing a robust template for correlating the ice-core records with other sediment records containing tephras (Lowe et al., 2008;Blockley et al., 2012;Abbott and Davies, 2012).In this work we focus on the identification of cryptotephra horizons (i.e.horizons that contain a low concentration of volcanic glass particles and are invisible to the naked eye).A sampling strategy was devised based on three criteria: -Ice spanning distinct sulfate peaks in the NGRIP CFA measurements following the approach of Davies et al. (2010) and Abbott et al. (2012).
-Intervals spanning rapid climatic transitions in the δ 18 O record.
-Intervals corresponding to age ranges of known tephra layers not yet identified in the ice cores.
In the sampled sections, ice strips with a cross section of approximately 2 cm 2 were cut from the 55 cm-long archive segments of the core stored at the University of Copenhagen.The sample preparation methodology followed that outlined in Davies et al. (2008Davies et al. ( , 2010) ) and Abbott et al. (2012), and involves centrifugation of the melted ice sample and preparation of the remaining particulate material onto ground glass slides.The samples were then embedded in epoxy resin and high-magnification light microscopy was undertaken to identify the presence of glass-shard particles.Samples found to contain volcanic material were then ground and polished using different grades of diamond suspension (9, 6 and 1 µm) to produce a polished surface suitable for analysing single shards by electron microprobe.Electron-probe microanalysis (EPMA) of the identified glass shards took place at the Tephra Analysis Unit based at the University of Edinburgh.A Cameca SX-100 electron microprobe with five vertical wavelength dispersive spectrometers was used to analyse oxide values for 10 major and minor elements within individual glass shards.The operating conditions followed those outlined by Hayward (2012).
Fifteen tephras common to the NGRIP and NEEM ice cores are identified and are outlined in Table 1.These are Icelandic in origin with 8 of tholeiitic basaltic composition and 7 of transitional alkali basaltic composition (Fig. 2).
Statistical comparisons testing the major-element correlations provide support for the matching of horizons between the cores (Table 1).Firstly, the statistical distance SD(D 2 ) test of Perkins (1995Perkins ( , 1998) ) which utilises the mean and standard deviations for the 10 elements that exceed 0.1 wt%, demonstrates that none of the correlations being tested are statistically different at the 99 % confidence interval as the D 2 values do not exceed the critical value of 23.21.Secondly, the use of the similarity coefficient (SC) function of Borchardt (1972), which compares the average values for the 7 major elements exceeding 1 wt%, highlights the strong similarities between the correlated deposits as SC values all exceed 0.95.This is usually taken as the threshold for defining an identical geochemical sample (Begét et al., 1992).

Methods and results
In the logging and analysis of ice cores, the core length is the fundamental parameter.The core is retrieved from a slightly inclined bore hole, and inaccuracies in the logged length accumulate along the core, especially in sections where drilling difficulties resulted in poor core quality; or in the brittle zone, where the core quality is often compromised by depressurisation of air bubbles within the ice.However, as all measurements on the core are made relative to the same logged depth scale, no such inaccuracies are of any relevance to the synchronization of the ice cores or timescale transfer.A potential offset between true depth and logged depth could be of relevance when records from the ice core are compared to in-hole measurements or modelled quantities.We have no direct way of estimating the magnitude of this offset, but the NEEM drilling and logging procedures are similar to those employed at NGRIP where the drilling and logging of two parallel ice cores down to a depth of 1372 m made it possible to evaluate the NGRIP2 logging precision to about 0.5 m at 1372 m depth (Hvidberg et al., 2002).As the estimated uncertainty on the logged depth therefore is small and varies slowly with depth, the effect on the comparison of true-depth ice flow modelling results and core-depth measurements is negligible.

Synchronization of the NEEM and NGRIP cores
To allow for the investigation of differences, including possible leads and lags, between the climate records of the two cores, the transfer of the timescale is based on reference horizons of non-climatic origin only.Ideally, only horizons that can be uniquely identified in both cores should be used.Extensive ongoing work on locating tephra in ice cores and geochemically characterizing these tephras is providing an ever-growing number of such horizons that tie different ice cores together with high confidence (Mortensen et al., 2005;Narcisi et al., 2005Narcisi et al., , 2012;;Coulter et al., 2012).However, the number of available tephra horizons is still far from sufficient to provide the sole basis for the transfer of GICC05 to the NEEM core.Instead, the synchronization is based on the identification of common patterns in the continuous and highly resolved ECM and DEP records with the tephra horizons providing an independent check of the synchronization and five additional match points.
Peaks in the ECM records originate from many sources with sulfuric acid of volcanic origin being a common source, while dips in the ECM record often are correlated to high NH + 4 concentrations, which are thought to be caused, for example, by episodes of high biomass burning activity (Fuhrer et al., 1996;Taylor et al., 1996).ECM peak heights vary significantly between cores probably due to differences in both transport and depositional conditions as well as variations in the background amount of alkaline dust and other impurities.Significant differences are seen even between replicate cores obtained from the same site (Wolff et al., 2005).However, even when individual peaks in themselves cannot be robustly matched between the cores, characteristic sequences of peaks and/or dips in the ECM signals can be used for robust climate-independent synchronization of ice cores (Rasmussen et al., 2008).In the section of NEEM below 1757 m depth, DEP data was used to support the matching.
The matching was performed by several investigators in parallel using the tool "Matchmaker" (available upon request from S. O. Rasmussen).All sections were matched independently by at least two investigators or groups of investigators.No significant differences between parallel matches were observed, although the number of match points differed markedly among investigators.In the second phase, a consensus set of match points was agreed upon by the investigators, excluding points only supported by one investigator.A total of 383 match points were identified between the NEEM and NGRIP1 cores with an additional 423 match points tying together the NEEM and NGRIP2 cores (of which 19 are in the zone of overlap between NGRIP1 and NGRIP2 measurements, leading to a total of 787 match points).The match points are provided in the data file accompanying this paper and a typical sequence of ECM and DEP data with match points is shown on Fig. 3.
During both phases of match-point identification, the Matchmaker tool is used to continuously evaluate whether the proposed match is glaciologically realistic (i.e.compatible with reasonable assumptions about the flow and accumulation regimes of the two core sites).As outlined in Rasmussen et al. (2008), the most useful diagnostic variable is the ratio r ab = (D a − D b ) / (d a − d b ) of annual layer thicknesses between neighbouring match points a and b found at depths D a and D b in the NEEM core and d a and d b in the NGRIP core.The values of r ab depend on both the relative accumulation rates of the cores and the amount of thinning induced by ice flow, and are expected to change only slowly (especially within a certain climatic state) as long as the spacing between a and b is large enough to eliminate the effect of short-term accumulation variability, surface redistribution, and uncertainties in the ECM depth scale (Andersen et al., 2006a).
The match points and the annual layer thickness ratios r ab are shown on Fig. 4 together with the location of the tephra horizons described in Sect.2.4.Long sections without match points occur during glacial stadial conditions where the relatively high dust levels make the ice alkaline, thus muting the ECM acidity record (Ruth et al., 2003).As seen in Fig. 4, the higher NEEM accumulation leads to thicker annual layers in NEEM than NGRIP down to a depth of about 700 m, after which the increased flow-induced thinning at NEEM leads to NGRIP annual layers being thicker by a factor of 1.5-2 during most of the glacial.
The (NGRIP depth, NEEM depth) relation has a small kink at the 2087.698m match point and clear kinks at the 2166.449m and 2195.630m match points.The first two kinks probably represent the boundaries between sections that have undergone different amounts of strain due to complex ice flow patterns that lead to overturned folds deeper in the ice (NEEM community members, 2013), but we believe that the ice stratigraphy is still undisturbed.Below the 2195.630m match point, the match between the cores is less robust, and there is a significant risk that the ice below this depth is stratigraphically disturbed.However, we do not have enough data of significant resolution at this point to shed light on this issue, and we therefore transfer GICC05modelext to NEEM down to a depth of 2203.597m (or 108.2 ka b2k) below which the ice core is known to be stratigraphically compromised, noting that the section below approx.2150 m (93.6 ka b2k) is less certain and that the section below the 2195.630m match point (107.0ka b2k) must be considered tentative.

The relation between tephra horizons and the ECM synchronization
As described in Sect.2.4 and Table 1, fifteen tephra horizons common to the NGRIP and NEEM cores have been located.All of these are fully consistent with the ECM-based synchronization.
The tephra particles are found in ice samples with a typical length of 15-20 cm, but their precise locations have not been determined by higher-resolution sampling and the tephra layers are not visible.In contrast, the ECM record has a much finer depth resolution and precision.More than half of the tephra horizons sit on top of (or very close to) one or more ECM-based match points, in which case only the more precise ECM-based match points are used for the timescale transfer.timescale transfer data set using the middle of the tephra sampling intervals (see Table 1 for details).

Transfer of the GICC05modelext timescale from NGRIP to NEEM
The transfer of the timescale is based on the ECM and DEP-based match points and the five tephra horizons located away from ECM-based match points.Using the match points' NGRIP depths and the NGRIP GICC05modelext (depth, age) relationship, the ages of the 787 ECM-based match points and 5 tephras are calculated.This set of mark-ers defines the (depth, age) relation for NEEM, named GICC05modelext-NEEM-1.The accuracy of the timescale at these points depends on three factors: -The NEEM timescale inherits the maximum counting error (MCE) of the NGRIP GICC05 timescale.
-Differences between the shape of peaks and inaccuracies in the depth registration of the ECM data set introduce synchronization uncertainty on the order of centimetres.The estimated synchronization uncertainty was estimated to 10 cm (1σ ) by Rasmussen et al. (2008), and here we tentatively estimate its magnitude by calculating the effect on the (NEEM depth, NGRIP depth) relation of removing every second match point.The results support the estimated synchronization uncertainty of 10 cm (1σ ), leading to timescale transfer uncertainties ranging from a few years to a maximum of a few decades at the deepest part of the record.
-Although we believe the set of match points to be robust, there is a risk that some sections have been erroneously matched up, leading to a larger systematic depth offset.This concern is particularly relevant for the NEEM ice below 2150 m (93.6 ka b2k).
As the MCE is typically 2 orders of magnitude larger than the matching uncertainty (when assuming no large systematic errors), we report GICC05modelext-NEEM-1 ages with the MCE uncertainty estimates only, but stress that observed phasing differences of up to a decade at the match-point depths could be artefacts from the timescale transfer.As seen on Fig. 4, there is a section across parts of the brittle zone (845-1026 m NEEM depth) with only few match points and a section of alkaline ice of GS-2 (1511-1595 m NEEM depth) with no match points at all.In addition, many shorter periods (mainly in the stadials) do not contain match points.Even in the absence of match points, we provide for convenience an interpolated (depth, age) relation for each 0.55 m depth increment, corresponding to the "bag" units used when cutting and packing the ice core.NEEM bag depths are first translated to NGRIP depths by linear interpolation using the (NEEM depth, NGRIP depth) relation established by match points and tephra horizons.By using linear interpolation, we assume a slowly varying ratio of annual layer thicknesses between the cores, which appears to be a reasonable assumption both from our understanding of accumulation variability and from the smoothness of the (NEEM depth, NGRIP depth) curve on Fig. 4. From the interpolated NGRIP depths, the ages are obtained from the NGRIP GICC05modelext (depth, age) relation.
Note that obtaining the bag ages directly from the 792 NEEM (depth, age) match points by interpolating the age linearly between match points gives a different result.This approach corresponds to assuming constant annual layer thicknesses between the NEEM match points, which, especially in sections with match points in interstadials only, is an unrealistic assumption given the large variations in accumulation rates across stadial-interstadial boundaries.

The precision of the interpolated timescale
As described in Sect.3.2, the NEEM and NGRIP ice cores are precisely aligned at the match points and the accuracy of the NEEM timescale at the match points is thus essentially quantified by the MCE of GICC05.However, the interpolation applied to provide a continuous NEEM timescale may introduce artificial offsets between the NEEM records and other records on the GICC05 timescale.These offsets are generally much smaller than the MCE and mainly matter when discussing the relative timing between events as observed in records from different cores.
As an estimate of the upper bound of the interpolation uncertainty, the effect of leaving out individual match points has been calculated.The interpolation was repeated using the entire set of match points except one (the leave-one-out approach).Going through all the match points in turn, interpolation differences of typically a few tens of centimetres or less are observed, with differences peaking at 95 cm when the excluded match point is the last point within GI-17 where the distance to the neighbouring match points is unusually large.
Using linear interpolation between match-point depths certainly is an approximation and it is not straightforward to evaluate its uncertainty.However, the smooth shape and small curvature of the (NEEM depth, NGRIP depth) curve seen on Fig. 4 illustrates that it is hard to justify more complex interpolation schemes.To test the influence of different interpolation schemes, interpolation was also performed using cubic spline interpolation.The change in interpolation method results in differences in the interpolated timescale of only 0-2 yr for the most recent 16 ka and no more than 20 yr for the period 96-69 and 60-23 ka b2k.These uncertainties are smaller than or equivalent to the synchronization uncertainty discussed in Sect.3.2.Exceptions are the long GS-2 match-point gap spanning 1511-1595 m NEEM depth and periods of similar duration without match points in GS-18 and GS-19.Here the possible offset could be larger, but we have little data from which to estimate this.Rasmussen et al. (2008) found unexpected offsets of a few metres in the (depth, depth) relations of the NGRIP, GRIP, and GISP2 cores across the GS-2 match-point gap interval, and although these may be related to ice divide migration or lee-side effects near the GRIP and GISP2 drill sites, it cannot be ruled out that similar issues exist for the NGRIP-NEEM synchronization.
Based on this and on visual inspection of the (NEEM depth, NGRIP depth) curve, we estimate that the likely maximum depth offset caused by interpolation across the GS-2 match-point gap is 1 m, corresponding to around 80 yr.However, the unusually strong accumulation variability across GS-2 could potentially introduce an interpolation-based offset significantly larger than this.Similarly, across the ∼ 5 kyr long GS-18 and GS-19 sections where the synchronization is based on only two match points within GI-18, an interpolation-based offset of 1 m (here corresponding to around 200 yr) is conceivable.

Estimating the thinning function
Past accumulation rates can be inferred from the observed annual layer thicknesses in an ice core if the strain history at the drilling site is known.accumulated vertical strain along the NEEM core, has been calculated using a Dansgaard-Johnsen (D-J) model (Dansgaard and Johnsen, 1969).
As discussed by NEEM community members (2013), the stratigraphy of the NEEM core has been compromised below 2203.6 m and, therefore, we do not attempt to model the strain in the deepest part of the core.A D-J model including basal melt and sliding (Buchardt, 2009) has been fitted to the (depth, age) relation of Sect.3.2.The basal melt and sliding have no straightforward physical interpretation given the disturbed and folded stratigraphy but has been included in the model in order to properly represent the flow regime of the NEEM area as reflected in the upper 2203.6 m of core.
The ice flow model is driven by a dynamical accumulation model, where the accumulation rate is calculated from the dated δ 18 O record through an exponential relation (Johnsen et al., 1995).The measured δ 18 O values have been corrected for changes in ocean δ 18 O (Bintanja and van de Wal, 2008), and an extended thickness history based on Vinther et al. (2009) is used.The inverse problem of determining the values of the unknown parameters in the ice-flow model and the accumulation model is solved using Monte Carlo analysis constrained by the observed (depth, age) relation.Though match points are identified back to 108 ka b2k, only depthage horizons younger than 85 ka b2k are used to constrain the flow parameters due to the accumulation ratio kinks discussed in Sect.3.2 and shown in Fig. 4, which cannot be adequately described using a D-J model.As the ice is disturbed and folded below 2203.6 m, the layers immediately above are likely to be influenced and adequate dating with the simple D-J model cannot be expected.
The Monte Carlo analysis reveals an average ice equivalent accumulation rate over the last 3 ka of 0.221 m yr −1 (1σ = 0.003 m yr −1 ), a typical glacial accumulation rate of 0.05 m in the stadials, and about twice this value in interstadials.However, it is not possible with this simple model to obtain a good correspondence between modelled and observed depth-age horizons in the Holocene.The modelled depths are too shallow for match points younger than ∼ 8 ka and too deep for Holocene points older than this.To further investigate this discrepancy, we look at the average annual layer thicknesses between the match points presented in Sect.3.1.These are shown in blue in Fig. 5 together with the annual layer thicknesses calculated for a constant Holocene accumulation rate of 0.221 m yr −1 using a thinning function calculated from the Monte Carlo estimated parameters (Fig. 5, red curve).A prominent kink is seen in the observed annual layer thicknesses around 1200 m depth corresponding to an age of approximately 7.9 ka b2k.This kink cannot be simulated with the simple D-J model, so in order to obtain a reasonable strain history for the ice above the kink, the Monte Carlo analysis is carried out constrained only by match points younger than 7.9 ka b2k.This leads to a better agreement between observed and modelled annual layer thicknesses for the period back to 7.9 ka b2k (Fig. 5, green curve).
The NGRIP and NEEM annual layer thicknesses show a parallel behaviour with early Holocene layers being relatively thin compared to what is expected from the model.Also, annual layer thicknesses from GISP2 seem to be somewhat thinner than expected from the model of Alley et al. (1993) between 11.7 and 9 ka b2k, whereas no such shift around 10-8 ka b2k is found in the Holocene annual layer thicknesses from GRIP (Vinther et al., 2006).
To further pursue the idea of a different early Holocene accumulation pattern, we compare observed mean annual layer thicknesses between match points with layer thicknesses derived by the fitted D-J model (black curve in Fig. 6).A large discrepancy is seen in the early Holocene, across GI-1, and during the late glacial.If the δ 18 O-inferred accumulation rates are reduced by 30 % during the glacial and the reduction in accumulation rate is linearly phased out from 30 % at 11.7 ka b2k to zero at 7.9 ka b2k, the match between the Monte Carlo-tuned D-J model and observed annual layer thickness improves dramatically (green curve in Fig. 6).
The origin of the changed accumulation-δ 18 O relationship in the early Holocene is unclear, and could be caused by changes in the accumulation-temperature relationship, changes in the temperature-δ 18 O relationship, limited ability of the flow model to represent a changing flow regime, or a combination hereof (Vinther et al., 2009).Possible Observed NEEM annual layer thicknesses (blue) and annual layer thicknesses inferred by applying the thinning function obtained from the D-J model to the accumulation rates calculated from the δ 18 O record (black).By reducing the calculated accumulation rates by 30 % in the glacial and phasing out the reduction linearly from 30 % at 11.7 ka b2k to zero at 7.9 kyr b2k, a significantly better fit to the observed annual layer thicknesses is obtained (green).
explanations include (i) ice ridge movements causing the need for the flow parameters in the D-J model to be adjusted, (ii) circulation and accumulation pattern changes caused by changes in Greenland ice sheet height and geometry and the gradually diminishing effects of the Laurentide ice sheet, and (iii) changes in δ 18 O caused by elevation changes and upstream effects.
The reconstructed accumulation history with the abovementioned forced accumulation reduction in the glacial and early Holocene is presented in Fig. 9.

Establishing the gas chronology by estimating age
Air bubbles are trapped at the firn-ice transition, leading to an ice age-gas age offset that is commonly referred to as age (Schwander and Stauffer, 1984).To obtain a gas chronology for the NEEM core, we model the evolution of age in the past using a coupled firn densification-heat diffusion model, with additional gas phase data to constrain the reconstruction (Fig. 7).We use δ 15 N of atmospheric nitrogen (N 2 ) where available as well as methane (CH 4 ) data obtained through a combination of discrete sampling (Mitchell et al., 2011;Rosen et al., 2013;Guillevic et al., 2013) and continuous flow analysis (Stowasser et al., 2012;Schupbach et al., 2009;Chappellaz et al., 2013).
The δ 15 N provides a strong constraint on both the timing and magnitude of abrupt temperature changes over Greenland through the imprint of thermal isotopic fractionation (Leuenberger et al., 1999;Severinghaus et al., 1998) and on past firn column thickness through the imprint of gravita-  (Chappellaz et al., 2013).tional isotopic fractionation (Sowers et al., 1992).CH 4 variations are in phase with Greenland temperature, with CH 4 lagging temperature by 0-70 yr (Huber et al., 2006;Landais et al., 2004;Severinghaus and Brook, 1999;Vallelonga et al., 2012).Thus, we obtain additional timing constraints by assuming the midpoint of the CH 4 transition to slightly lag the midpoint in the δ 18 O transition at the abrupt onset of Greenland interstadials (GI); for simplicity we use a constant lag of 30 yr for most transitions, with the exception of interstadials 9-11 where available δ 15 N data indicate no lag.In several cases abrupt CH 4 and δ 18 O features are observed within an interstadial, providing additional timing constraints.The CH 4 tie points for interstadials 5-11 are indicated by white circles in Fig. 7.
For the gas chronology presented here, δ 15 N data were available for the last glacial termination (∼ 15-11 ka b2k) and for interstadials 8-11 and 21; reliable CH 4 constraints were available for the glacial termination and interstadials 2-22.The depth constraints from the gas data were converted to age constraints through the NEEM ice chronology as presented in Sect.3.2.A final constraint is the modern age of 188 yr as derived from firn air measurements (Buizert et al., 2012).
Due to the many modelling uncertainties, the data-derived timing constraints must be considered more accurate than the modelled age.Because of the strong variability in

S. O. Rasmussen et al.: A first chronology for the NEEM ice core
Greenland temperature and accumulation during the glacial period, age can change by hundreds of years in between the tie points, precluding the use of simple interpolation schemes.The philosophy behind our age reconstruction is to use densification models to obtain a realistic, physically meaningful interpolation between the data-derived tie points.Past changes in accumulation rate and surface temperature are the main drivers of age variability, but both input functions are poorly known back in time.We adjust the accumulation (A) and temperature (T ) input to minimise the RMS difference between the modelled ages and the timing constraints.
We use two different firn densification models in this study, with supporting information given in the Appendix.The first is the model by Goujon et al. (2003), which is a dynamical version of the 1-D densification model by Arnaud et al. (2000).The heat diffusion in the ice is calculated across the entire ice depth using a simplified version of the model from Ritz (1989).The model uses a time step of one year for both densification and heat diffusion calculation.The spatial resolution in the top 150 m is 0.25 m; it decreases with depth to 25 m at bedrock.The lock-in depth (LID), below which the air is no longer exchanging with overlaying layers, is defined as the depth where the ratio of closed to total porosity reaches 0.13.This LID definition based on the porosity of the ice is in agreement with the LID density threshold proposed by Schwander et al. (1997).A constant convective zone of 2 m is used.This model is validated for present-day conditions in Greenland and Antarctica (Landais et al., 2006;Arnaud et al., 2000).
The second model is a dynamical version of the Herron-Langway model (D-HL) (Herron and Langway, 1980).The model has 0.5 m spatial resolution along the depth axis, with the model domain extending down to 1000 m.The heat diffusion model uses implicit Crank-Nicolson time stepping, a zero gradient boundary condition at 1000 m depth and an advective heat flux based on the 1-D ice flow model presented in Sect.3.3.The densification and heat diffusion models use 1.0 and 0.2 yr time steps, respectively, with isochrones tracked downwards every fifth year.Lock-in density is set at 14 kg m −3 below the mean close-off density from Martinerie et al. (1994) as suggested by Schwander et al. (1997); lock-in gas ages are calculated using the parameterization of Buizert et al. (2013).Following recent work by Hörhold et al. (2012), the model includes an empirical softening of ice that scales with the logarithm of the Ca 2+ concentration, where we use Ca 2+ data from Greenland summit (Fuhrer et al., 1993) as a measure of the dust content.The softening is achieved by multiplying the thermal activation energy by the factor 1 − αln β × [Ca 2+ ] (J.Freitag, personal communication, 2012); the fitting factors α = 0.0021 and β = 3 ppb −1 were calibrated for modern NEEM firn where accurate T , A, age and [Ca 2+ ] data are available.Without including this empirical dust softening term, the densification model does not re-produce the present-day age of 188 yr correctly (Buizert et al., 2012).While yielding a more accurate present-day age, including the dust effect does not influence the results in the deep part where age is constrained by δ 15 N and CH 4 data.In Appendix A we model age for the last deglaciation (20-10 ka b2k) using α = 0, α = 0.0021, and α = 0.0042, and find that the results are identical within the uncertainty (Fig. A2).
The initial estimate for T (t) is based on the δ 18 O palaeothermometer, using a sensitivity of 0.562 ‰ K −1 and a sea-water δ 18 O correction (Waelbroeck et al., 2002).The isotope sensitivity used here is found by δ 15 N calibration of the δ 18 O palaeo-thermometer over interstadials 8-10 at NEEM (Guillevic et al., 2013).The initial estimate for the accumulation is based on a combination of straincorrected annual layer thicknesses (similar to the red curve of Fig. 9), and a relationship with δ 18 O for the deeper core where the strain correction becomes unreliable: A = exp 0.144 × δ 18 O + 3.245 .These initial A and T reconstructions are subsequently adjusted to optimize the fit to the data; this is done manually for the Goujon model, and with an automated gradient method in the D-HL model.Adjusted A and T profiles, as well as modelled age values are shown for both models in Appendix A, while the final accumulation histories for the D-HL and the Goujon models are shown on Fig. 9.The D-HL model provides a better fit to both the δ 15 N and CH 4 age constraints; we attribute this to the automated calibration routine, rather than differences in the model physics.The root-mean-square (RMS) offset between model and CH 4 tie points is 19 yr (66 yr) for the D-HL (Goujon) model, while the average tie point uncertainty is 63 yr.The RMS offset to the δ 15 N data is 0.011 ‰ (0.018 ‰) for the D-HL (Goujon) model, compared to an analytical precision on the order of 0.007 ‰.Consequently we use the D-HL model result for periods where we have sufficient timing constraints; while during two periods without constraints (108-90 and 23-17 ka b2k), output from both models is averaged.The combined age curve is shown in Fig. 8 together with the age constraints from CH 4 .The age uncertainty σ age is estimated in two ways.First, at the CH 4 tie points we set σ age equal to the uncertainty in that tie point (the error bars shown on Fig. 8) plus the absolute value of the modeltie point mismatch.In between the tie points we set σ age to the linearly interpolated uncertainty of the two adjacent tie points plus 0.05 × t, with t being the distance to the nearest tie point in years; after the last tie point at 88 ka b2k we keep σ age constant.As a second estimate we use the absolute value of the Goujon minus D-HL model difference.For the glacial period (108-11.7 ka b2k) we use the larger of the above two estimates as the final age uncertainty (shaded area on Fig. 8).In the Holocene we tentatively set σ age at 9 ka b2k to 30 yr, and again use linear interpolation to find σ age at other times.

Reconstructing past accumulation rates
From the ice flow and firn modelling, we obtain four different estimates of past accumulation rates: -A δ 18 O-based accumulation profile entering the ice flow model of Sect.3.3 is determined from the Monte Carlo analysis.
-Observed annual layer thicknesses (here we use 20 yr means) from NGRIP can be transferred to NEEM depths using the match points as described in Sect.3.2, and by correcting for the ice-flow-induced thinning using the modelled strain (Sect.3.3), past accumulation rates are inferred.
-As outlined in Sect.3.4, past accumulation rates are estimated by two different firn models by adjusting accumulation and temperature histories to optimize the fit to the age constraints.
The four different estimates are shown back to 85 ka b2k on Fig. 9, and it is clear that there is considerable disagreement between the results especially in the 55-30 ka b2k interval, with the estimates based on firn modelling being typically 30 % lower than the estimates based on ice flow modelling.We want to stress that without knowledge of past convective zone thicknesses and surface densities, and in the absence of high resolution δ 15 N data for most of the core, firn modelling-derived past accumulation rates are fundamentally under-constrained.Similarly, the ice flow modelling is based on a simple one-dimensional D-J model at a location where the basal ice is folded and disturbed.However, we believe that the combined results provide a realistic range of past NEEM accumulation rates, and we therefore provide the mean and 2σ envelope of the estimates (black curve and grey shaded area in Fig. 9).

Summary and perspectives
The GICC05modelext timescale has been applied to the NEEM ice core by transfer from NGRIP using 787 horizons identified in the ECM and DEP records.Tephra layers confirm this synchronization and add five additional match points for the timescale transfer.The gas record has been dated by determining age from firn modelling with constraints from δ 15 N, CH 4 and modern-day age values, while the flow-induced strain history is determined from an ice flow model.The accumulation reconstructions from the ice flow and firn models show significant differences, in particular during the stadials in the middle part of the glacial (55-30 ka b2k).The firn model reconstruction can be improved when δ 15 N data becomes available, while the ice flow model based results can be improved with improved ice flow models and additional knowledge on the evolution of the ice sheet.In addition, steps should be taken to treat ice flow and firn processes in an integrated way with the use of common accumulation and temperature input profiles.
The presented GICC05modelext-NEEM-1 timescale is intended as a first chronology for the NEEM ice core facilitating the analysis of existing and forthcoming data sets from the NEEM core.When more tephra horizons and highresolution CFA data from NEEM become available, a more detailed synchronization will be possible, hopefully also providing more match points within stadial periods and thereby increasing the precision of the timescale transfer.Future investigations will also show if NEEM data can be used to refine the annual layer identification in ice cores from Greenland, e.g. in the period 10.3-7.9 ka b2k.

Data access
The following data sets are made available at www. iceandclimate.dk/data and WDC Palaeo: -ECM data sets from NGRIP and NEEM in 1 cm resolution.
-NEEM DEP data from 1757 m depth downwards and NGRIP DEP data.
-The match points (including the five tephra horizons) used for the timescale transfer.
-The NEEM (depth, age) relation in 0.55 m ("bag") resolution for ice and gas.
-Modelled total strain and accumulation in 0.55 m resolution.

Appendix A
Figure A1 shows the temperature T and accumulation rates A used in the age modelling.The densification models are essentially run as inverse models, where one is looking for the A and T that optimize the fit to the CH 4 tie points and δ 15 N data.The models start from an initial guess of A and T , given by the grey curves.A description of the initial curves is given in Sect.3.3.Modifications to the A and T curves are done manually for the Goujon model and using an automated gradient method for the D-HL model.For most of the core, only CH 4 tie points are available, which leaves the A and T under-constrained.Therefore no firm conclusions should be drawn from these reconstructions.
Figure A2 illustrates the influence of changing the parameter α quantifying the importance of the dust softening term described in Sect.3.4.

Fig. 2 .
Fig. 2. Biplots showing the tephra horizons matched between the NEEM and NGRIP cores.The tephras are numbered according to Table 1.Correlative tephra horizons are represented by the same shape and colour symbols, with NEEM layers represented by open symbols and NGRIP layers by filled symbols.(A) biplots for the tephra layers included in the timescale transfer.(Ai) plot of SiO 2 vs. CaO and (Aii) plot of FeO/MgO vs. TiO 2 .(B) biplot of CaO vs. TiO 2 for the remaining tholeiitic basalt layers and (C) biplot of CaO vs. TiO 2 for the remaining transitional alkali basalt layers.

Fig. 3 .
Fig. 3. Example of ECM and DEP data across a synchronized section of the NEEM and NGRIP cores.The grey bands mark the match points used for the timescale transfer.

Fig. 4 .
Fig. 4. ECM/DEP-based match points between the NEEM and NGRIP1 (red) and NGRIP2 (blue) ice cores and the tephra horizons (black squares).In green (right axis), the NGRIP/NEEM annual layer thickness ratio r ab is shown.The insert shows the details of the match points and r ab over the lowest 150 m of the core, where kinks in the (NEEM depth, NGRIP depth) relation and corresponding jumps in r ab most likely indicate areas subject to differential strain.

Fig. 5 .
Fig.5.Mean NEEM annual layer thicknesses between the match points presented in this work (blue dots) and predicted Holocene annual layer thicknesses calculated from a Dansgaard-Johnsen model with constant accumulation rate of 0.221 m yr −1 and flow parameters determined from Monte Carlo analysis constrained by observed depth-age horizons back to 85 ka b2k (red line).The green line shows the predicted annual layer thicknesses from a D-J model with flow parameter values determined from Monte Carlo analysis constrained only by depth, age horizons younger than 7.9 ka b2k.The transition from the glacial to the Holocene (11.7 ka b2k) is marked with a horizontal grey line.High-frequency variability in the observed annual layer thicknesses caused by closely spaced match points in the top part of the core has been removed.
Fig.6.Observed NEEM annual layer thicknesses (blue) and annual layer thicknesses inferred by applying the thinning function obtained from the D-J model to the accumulation rates calculated from the δ 18 O record (black).By reducing the calculated accumulation rates by 30 % in the glacial and phasing out the reduction linearly from 30 % at 11.7 ka b2k to zero at 7.9 kyr b2k, a significantly better fit to the observed annual layer thicknesses is obtained (green).

Fig. 7 .
Fig. 7. NEEM data used in age reconstruction at the onset of interstadials 8-11.(a) δ 18 O of precipitation as a proxy for site temperature with 0.55 m averages in light blue, and a 5-point running mean in dark blue; (b) δ 15 N of N 2 from Guillevic (2013) (black dots) and model output (green curve); (c) methane data from continuous flow analysis: PICARRO instrument (orange) and SARA instrument (red)(Chappellaz et al., 2013).

Fig. 8 .
Fig. 8. Evolution of NEEM age with timing constraints from CH 4 data (green curve and shaded area indicating the estimated uncertainty).Also shown are the results from the Goujon (blue) and dynamical Herron-Langway (orange) models and age constraints derived from CH 4 data (black dots).

Fig. 9 .
Fig. 9. Accumulation reconstructions from the two firn models (Goujon in blue and dynamical Herron-Langway in orange), from the δ 18 O-based parameterization with 30 % reduction in glacial accumulation that is used in the ice flow model (red), and inferred by transferring 20 yr mean annual layer thicknesses from NGRIP to NEEM and taking modelled strain into account (green).The black curve and grey envelope show the mean and 2σ of the four individual reconstructions.
Fig. A1.Temperature (T , upper) and accumulation (A, lower) input into the densification models: results from the Goujon (blue) and dynamical Herron-Langway (orange) models are shown together with the grey curves, which show the initial T and A estimates as described in the text.

Table 1 .
Bourne et al. (2013)05) horizons identified in the NEEM and NGRIP ice cores.Horizons used for the timescale transfer are marked by a , while the other horizons are consistent with the synchronization but not used for the timescale transfer as they are located too close to better-resolved ECM match points.Previously published major element data for NGRIP horizons 1-3 are fromMortensen et al. (2005)and horizons 5 and 7 are fromDavies et al. (2010).Major element results for horizons 11-15 are reported inBourne et al. (2013)and the remaining data sets will be available in forthcoming publications.bThetephraatNGRIP2066.95m was originally identified byDavies et al. (2010)and new geochemical data reported byBourne et al. (2013)and used in this work demonstrate that this horizon correlates to NEEM 1756.90-1757.10m.cMajorelementresults for NEEM 1764.05-1764.25 m are indistinguishable from two closely-timed tephras inNGRIP (2077NGRIP ( .90-2078NGRIP ( .01 m and 2078NGRIP ( .30-2078.37 m) (see Fig.2) and based solely on the geochemical signatures it is uncertain which one is the NGRIP correlative.However, only the correlation to NGRIP 2077.90-2078.01m is consistent with the ECM-based match.The columns SD(D 2 ) and SC provide the statistical distance and similarity coefficient values, respectively (see Sect. 2.4 for further details).