How old is the Ordovician–Silurian boundary at Dob’s Linn, Scotland? Integrating LA-ICP-MS and CA-ID-TIMS U-Pb zircon dates

Abstract

Sedimentary rocks exposed at Dob’s Linn, Scotland, have significantly influenced our understanding of how life evolved over the Ordovician to Early Silurian. The current interpreted chronostratigraphic boundary between the Ordovician and Silurian periods is a Global Boundary Stratotype Section and Point (GSSP), calibrated to 443.8 ± 1.5 Ma (Hirnatian–Rhuddanian age), based on biostratigraphic markers, radioisotopic dates and statistical modelling. However, challenges arise due to tectonic disturbances, complex correlation issues and the lack of systematic dating in Ordovician–Silurian stratigraphic sections. Here, hundreds of zircon grains from three metabentonite ash horizons were dated using Laser Ablation Inductively Coupled Plasma Mass Spectrometry (LA-ICP-MS). A subset of the grains were re-analyzed using Chemical Abrasion Isotope Dilution Thermal Ionization Mass Spectrometry (CA-ID-TIMS). We present a high-precision CA-ID-TIMS ²³⁸U-²⁰⁶Pb weighted mean date of 440.44 ± 0.55/0.56/0.72 Ma (±analytical/with tracer/with U-decay constant) for the Coronagraptus cyphus biozone. However, the study reports younger, and in certain cases, older LA-ICP-MS zircon dates within the Coronagraptus cyphus, Akidograptus ascensus and Dicellograptus anceps zones, suspected as being influenced by Pb loss and LA-ICP-MS matrix mismatch. The study reports concerns about the suitability of Dob’s Linn as a GSSP section and examines various LA-ICP-MS maximum depositional age (MDA) approaches, suggesting the use of the TuffZirc date and the youngest mode weighted mean (YMWM) as suitable MDA calculations consistent with CA-ID-TIMS results.

1. Introduction

Evaluating the boundary between the Ordovician–Silurian periods and interpreting the timing and duration of environmental and biological changes requires precise and accurate dating of stratigraphic sections that contain rocks of these time frames. This specific boundary is fundamental to comprehending how life appeared and radiated on this planet as well as understanding the timing of the Late Ordovician Mass Extinction (LOME) (Lenton et al. 2012; Wallace et al. 2017; Servais et al. 2019; Dahl et al. 2021). Calibrating relative timescales with isotopic dating of igneous rocks has been an ongoing task since the early days of radiometric dating (Holmes, 1911). Biostratigraphic boundary ages are continually revised with new methods, concepts and studies (Mattinson, 2013; Gradstein & Ogg, 2020). The development of accurate and precise zircon U-Pb dating methods has revolutionised the calibration of many parts of the geologic timescale (Bowring et al. 2006; Schoene et al. 2013; Spencer et al. 2016). The Ordovician and Silurian, however, suffer from a lack of systematic dating of volcanic lavas, breccias and ashes interstratified with biostratigraphically dated sediments. Furthermore, many local biostratigraphic schemes for different areas cannot be accurately correlated between marine and non-marine sections. Thus, the North American, British and Scandinavian schemes suffer from a number of correlation problems, and the Mediterranean and North Gondwanan schemes, and it is complicated to relate to the standard Series and Stages (Sweet & Bergström, 1984; Berry, 1987; Finney, 2005; Fortey, 2011).

Laser Ablation Inductively Coupled Plasma Mass Spectrometry (LA-ICP-MS) and Secondary Ion Mass Spectrometry (SIMS) allow for rapid U-Pb dating of zircons to determine provenance and maximum depositional ages (MDAs) (Table 1). Some studies show that LA-ICP-MS and SIMS methods have systematic biases in ²³⁸U-²⁰⁶Pb zircon dates relative to those obtained using the more precise yet destructive Chemical Abrasion Isotope Dilution Thermal Ionization Mass Spectrometry method (CA-ID-TIMS) (Mattinson, 2005; Allen & Campbell, 2012; Crowley et al. 2014; Marillo-Sialer et al. 2014; Von Quadt et al. 2014; Watts et al., 2016; Catlos et al. 2021). For these reasons, understanding the U-Pb zircon system and potential systematic biases is crucial to produce accurate dates for MDAs of sections of interest.

Table 1. List of commonly used maximum depositional age (MDA) methods modified from Sharman and Malkowski (2020)

Here, we report a high-resolution zircon ²³⁸U-²⁰⁶Pb CA-ID-TIMS date for the Coronagraptus cyphus biozone and test various LA-ICP-MS MDA calculations to determine suitable MDA approaches for the Coronagraptus cyphus, Akidograptus ascensus and Dicellograptus anceps zones from the tectonically disturbed metabentonites encompassing the Ordovician–Silurian boundary at Dob’s Linn, Scotland (Fig. 1). We compare CA-ID-TIMS and LA-ICP-MS zircon U-Pb dates from three graptolite zones showing significantly younger and, in some cases, older LA-ICP-MS zircon dates than anticipated. The discrepancies between CA-ID-TIMS and LA-ICP-MS with the younger and older LA-ICP-MS zircon dates indicate that such dates need careful evaluation and can be variously interpreted as a result of Pb loss, matrix mismatch and/or potential biostratigraphic misplacement, bringing doubt to the validity of Dob’s Linn as a Global Boundary Stratotype Section and Point (GSSP) reference section.

Figure 1. (a) Generalized map of the United Kingdom showing Dob’s Linn study area (red) in the Southern Scottish Uplands. (b) Geological timescale with respective graptolite zones. The red line presents the base of the Akidograptus ascensus biozone representing the Ordovician–Silurian boundary.

2. Geochronology Background

Zircon is a mineral that is resistant to weathering and thus often used to date the MDA of sedimentary sections (Carroll, 1953; Balan et al. 2001; Finch & Hanchar, 2003). U-Pb zircon geochronology is considered the optimal radioisotopic dating approach because two decay schemes generate two independent chronometers that can be cross-validated over geologic time. The two independent radioactive decay schemes consist of ²³⁵U-²⁰⁷Pb and ²³⁸U-²⁰⁶Pb, each with a different half-life, permitting identification of inherited domains and open-system behaviour (i.e., Pb loss) (Bowring et al. 2006; Corfu, 2013; Schoene, 2014). Three U-Pb dating methods can be used to date lavas and ashes; however, accuracy and precision vary significantly depending on the dating technique (Condon & Bowring, 2011; Spencer et al. 2016). As seen in Fig. 2, these approaches sample different portions of the zircon and yield different ranges of precision (Bowring et al. 2006; Condon & Bowring, 2011; Spencer et al. 2016). LA-ICP-MS is a high-speed and cost-effective dating technique with moderate precision, but depending on preparation methods, such as whether the unknown and standard zircons are annealed or not, can produce 2σ analytical precision of 1–8% (Von Quadt et al. 2014; Schaltegger et al. 2015; Ver Hoeve et al. 2018). Zircon grains are analyzed with a 10–60 µm spot size and 5–20 µm laser depth at a rate of 20 second–4 minutes per analysis (Bowring et al. 2006; Mako et al. 2021). SIMS, which includes Sensitive High-Resolution Ion Micro Probe (SHRIMP), is a rapid technique with a 2σ precision of 1–5%, 10–20 µm spot size, <2 µm analysis depth and a rate of 10–30 minutes per analysis (Bowring et al. 2006; Schaltegger et al. 2015; Tichomirowa et al. 2019). The most precise and accurate technique is Isotope Dilution Thermal Ionization Mass Spectrometry (ID-TIMS), with an additional chemical abrasion (CA-ID-TIMS) option capable of removing radiation-damaged and Pb loss domains in zircon grains. ID-TIMS requires several days of preparation in clean chemistry lab environment and takes 5–6 hours per mass spectrometric analysis with a 2σ age precision of ≤ 0.3%, while CA-ID-TIMS improves the accuracy of the dates by eliminating the effects of Pb loss and produces 2σ age precisions of ≤ 0.1% (Mattinson, 2005; Bowring et al. 2006; Schaltegger et al. 2015). Isotope dilution using a well-calibrated isotopic tracer eliminates the dependence on standard measurements and potential matrix mismatches that limit the accuracy and precision of spot analyses (LA-ICP-MS and SIMS; Bowring et al. 2006).

Figure 2. Generalized comparison of zircon U-Pb dating methods. Modified from Bowring et al. (2006).

Although accuracy and precision from LA-ICP-MS, SIMS and CA-ID-TIMS vary, each of these dating techniques has their advantages and limitations for establishing sedimentary MDAs. CA-ID-TIMS, a time-consuming, costly and destructive technique, functions best with individually dated zircons when accuracy and the highest precision are required. These evaluations pair well with cathodoluminescence (CL) imaging and pre-screening with LA-ICP-MS or SIMS to target the youngest autocrystic grain population from a temporally distinct magmatic pulse and prevent the inclusion of antecrysts formed from an earlier magma pulse, or xenocrysts included from older host rock during younger magmatic pulses (Rossignol et al. 2019; Zellmer, 2021). Additionally, the number of zircon grains per sample from the youngest age mode is significant as there is no assurance that these grains will endure the destructive chemical abrasion process and fully dissolve altogether with the radiation damaged Pb loss zones. Rapid, cost-effective analytical techniques such as LA-ICP-MS and SIMS using laser-coupled plasma or ion bombardment possess a high spatial resolution ideal to target specific domains in samples with abundant quantities of zircon grains. Zircon grains are polished and CL imaged to avoid inherited cores or potential metamict zones. Alternatively, zircons can be depth profiled, providing core-rim spatial information and spread in uranium concentrations. Since zircons are not polished and no CL images are acquired, the depth profile method is not optimal for complex grains with abundant growth zone history (Marsh & Stockli, 2015; Rasmussen et al. 2019). However, systematic biases with LA-ICP-MS and SIMS, such as Pb loss and the matrix effect between unknown and standards, are suggested to be the driving mechanism producing discrepancies across dating techniques (Bowring & Schmitz, 2003; Andersen et al. 2019).

In the last decade, LA-ICP-MS studies indicate systematic biases with ²³⁸U-²⁰⁶Pb zircon dates relative to CA-ID-TIMS and within the LA-ICP-MS technique itself, varying between laboratories. Different laboratories present consistently young or older dates for the identical sample using the same calibration standards as a result of the matrix effect (Marillo-Sialer et al. 2014). Matrix mismatch is a recognized systematic bias with LA-ICP-MS U-Pb zircon geochronology that is yet to be entirely comprehended. The three primary factors associated with the matrix effect originate with the zircon grains (unknowns), zircon standards and mass spectrometer ablation conditions (Jackson et al. 2004; Allen & Campbell, 2012; Marillo-Sialer et al. 2014). It is difficult to obtain identical behaviour between unknown and standard zircons for various reasons, including differences in grain sizes, radiation damage (alpha dose) and maintaining equal laser beam conditions, including spot size, focus, ablation rates and integration time (Jackson et al. 2004). Von Quadt et al. (2014) suggest that the physical condition of the unknown zircon grains and the utilized standards are the underlying cause of downhole fractionation of Pb from U, resulting in the matrix effect. According to Marillo-Sialer et al. (2014), the primary limitation of LA-ICP-MS is the requirement of the same behaviour between standards and unknowns during analysis. Allen and Campbell (2012) propose that the mechanism driving the matrix effect is the difference between the alpha dose between unknown and standard zircons, thus generating LA-ICP-MS ²³⁸U-²⁰⁶Pb zircon dates younger or older relative to ID-TIMS due to fractionation.

Additionally, factors such as tectonics and hydrothermal alteration can increase radiation damage accumulation experienced by zircon grains, consequently producing metamict Pb loss domains on a case-by-case basis depending on the geologic history and location of the unknown zircon grains (Schoene, 2014). Because radiation damage in zircon grains can vary extensively, it is challenging to utilize a well-characterized zircon standard identical to any possible unknown zircon grains (Jackson et al. 2004). However, the concerns of matrix mismatch induced by zircons affected by alpha decay radiation damage domains and spontaneous fission in the crystal lattice can be minimized by annealing both zircon standards and unknowns prior to spot analysis (Allen & Campbell, 2012; Solari et al. 2015; Ver Hoeve et al. 2018). For example, the well-characterized zircon standards GJ-1 and Plesovice contain metamict sectors that if annealed can improve LA-ICP-MS accuracy and precision (Jackson et al. 2004; Sláma et al., 2008; Frei & Gerdes, 2009). According to Ver Hoeve et al. (2018), LA-ICP-MS downhole fractionation is one of the principal setbacks in optimizing both precision and accuracy. While thermally annealing unknown grains improves accuracy, if the standards are also annealed then precision can be improved by minimizing downhole fractionation and matrix mismatch by obtaining as close as possible identical behaviour between unknown and standard zircons (Ver Hoeve et al. 2018).

3. Geological Background

The Dob’s Linn locality is one of several tectonically disturbed sections of the Moffat Shale Group of Southern Scotland with many intricate minor- and large-scale faults, isoclinal folding, and unresolvable thinning and thickening of strata (Williams, 1988). As shown in Fig. 3, the Dob’s Linn Ordovician–Silurian boundary GSSP outcrop’s bedding is positioned in a vertical direction due to the severe tectonism in southern Scotland. The widespread tectonic activity is associated with the forming of the Caledonian mountains that started during the Early Ordovician (475 Ma) and fully formed by the Late Silurian (425 Ma) (Fig. 4) (Chew & Strachan, 2014). The Moffat shale consists of a pelagic mudrock succession deposited in oceanic, forearc or back-arc environments in Late Ordovician to Early Silurian times (Fig. 4) (Morris, 1987; Stone et al. 1987; Merriman & Roberts, 1990).

Figure 3. Moffat Shale (Birkhill Shale member) at Dob’s Linn’s GSSP Linn trench outcrop showing vertical stratigraphy, graptolite horizons and metabentonite horizons that are deformed as a result of extensive tectonic activity. The red line displays the Ordovician–Silurian boundary with an accepted age of 443.8 ± 1.5 Ma (Cohen et al. 2013; Cohen et al. 2022). The orange dashed line shows sample 19DL09 in the study, a metabentonite horizon in the Akidograptus ascensus zone.

Figure 4. Paleogeographic reconstruction during the Late Ordovician (443 Ma). Closing of the Iapetus Ocean forming volcanic arcs (fore-arc and back-arc) near subduction margins of Laurentia, Baltica and Avalonia, producing widespread tectonic activity and forming the Caledonian mountains (after Huff et al. 2010; Chew & Strachan, 2014).

The Upper Ordovician–Lower Silurian sediments consist of the 48-metres-thick Hartfell Shale Formation subdivided into the Lower and Upper units containing scarce metabentonite horizons. The Lower Hartfell Shale is composed of primarily black mudstone coarsening upwards to cherty and silty mudstone. The Upper Harfell Shale is characterized by laminated and bioturbated grey mudstone (Williams, 1983; Batchelor & Weir, 1988). The Birkhill Shale Formation is 45 metres thick and likewise sectioned into lower and upper units, with continuous metabentonite successions and a sharp contact between the units. The Lower Birkhill Shale is a black mudstone transitioning to brittle cherty mudstone and with blocky morphology. The Upper Birkhill Shale is a black mudstone transitioning to grey–green mudstone (Batchelor & Weir, 1988; Merriman & Roberts, 1990).

Dob’s Linn is considered by some as a significant location due to the appearance of critical graptolite transitions during the Ordovician (485.4–443.8 Ma) and Silurian (443.8–419.2 Ma) periods in the Moffat Shale Group (Fig. 3) (Cocks, 1985, 1988; Williams, 1988; Gradstein et al. 2020). This chronostratigraphic boundary was first dated based on biostratigraphic distributions of graptolites (Fig. 1b) (Carruthers, 1858; Nicholson, 1867; Lapworth, 1878). These small, aquatic colonial animals are “unrivaled in the Early Palaeozoic” in terms of subdividing relative time (Zalasiewicz, 2001: 240) and thus are widely used to correlate sedimentary sections that contain them regionally and globally (Koren’ & Rickards, 1979; Williams, 1983). However, the use of correlating these early organisms has been problematic due to local evolutionary provincialism and convoluted age interpretations with the correlation of fluvial and marine deposits with a standard geologic timescale (Berry, 1987; Finney & Chen, 1990; Pogson, 2009; Brookfield et al. 2021). Additionally, the Ordovician–Silurian boundary age at Dob’s Linn has been estimated by several radioisotopic dates with varying precisions, calculated by spline fitting interpolation from units stratigraphically above and below the boundary (Tucker et al. 1990; Hu et al. 2008; Schmitz & Ogg, 2020).

4. Biozone ages

4.a. Coronagraptus cyphus Biozone

The Coronagraptus cyphus biozone located in the Lower Birkhill Shale Formation constrains the end of the Early Silurian Rhuddanian Stage (Fig. 1b) (Ross et al. 1982; Tucker et al. 1990; Gradstein et al. 2020). A zircon fission-track date of 437 ± 10 Ma was initially reported for the Coronagraptus cyphus zone in Dob’s Linn (Ross et al. 1982). Elsewhere, hornblende from the Coronagraptus cyphus zone in the Descon Formation in Esquibel Island, Alaska, produced a ⁴⁰Ar-³⁹Ar date of 442.6 ± 5.0 Ma (Lanphere et al. 1977; Ross et al. 1982; Kunk et al. 1985; Schmitz & Ogg, 2012). Utilizing mechanically (air) abraded zircon ²³⁸U-²⁰⁶Pb ID-TIMS dates, Tucker et al. (1990) produced a biozone age of 439.57 ± 1.33 Ma from the weighted mean of six multigrain zircon fractions from Dob’s Linn, Scotland (Table 2) (Schmitz & Ogg, 2020).

Table 2. Compilation of radioisotopic dates and statistical approaches from previous studies that estimate the ages for graptolite biozones at or near the Ordovician–Silurian boundary

References 1. Ross (1984) 2. Tucker et al. (1990); 3. Schmitz & Ogg (2012); 4. Gradstein et al. (2012); 5. Ling et al. (2019); 6. Cohen et al. (2022); 7. Gradstein et al. (2020); 8. Ross (Ross et al. 1982); 9. Lanphere et al. (1977); 10. Kunk et al. (1985).

4.b. Akidograptus ascensus Biozone

The Akidograptus ascensus biozone located in the Lower Birkhill Shale Formation has been interpreted to define the Ordovician–Silurian boundary at Dob’s Linn, Scotland (Fig. 1b) (Rong et al. 2008; Gradstein et al. 2020). The Akidograptus ascensus zone is dated using various data points including graptolites and stratigraphically upper and lower zircon ID-TIMS U-Pb dates with spline fitting interpolation generating random replications with the input data and validated with a smoothing factor value producing a straight-line fit (Agterberg et al. 2020). However, the current age interpretations of 443.8 ± 1.5 Ma or 443.1 ± 0.9 Ma are calculated by spline fitting interpolation from U-Pb multigrain zircon fractions stratigraphically above (Coronograptus cyphus zone) and below (Dicellograptus anceps zone) the Akidograptus ascensus zone. Gradstein et al. (2020) interpolated the Ordovician–Silurian boundary age to 443.1 ± 0.9 Ma from Katian and Rhuddanian ID-TIMS ²³⁸U-²⁰⁶Pb zircon dates from Dob’s Linn and a Hirnantian SHRIMP ²³⁸U-²⁰⁶Pb date from South China (Tucker et al. 1990; Hu et al. 2008; Ogg et al. 2016; Schmitz & Ogg, 2020). The International Commission of Stratigraphy interpolated the age of the boundary to 443.8 ± 1.5 Ma only using Dob’s Linn’s Katian and Rhuddanian ID-TIMS ²³⁸U-²⁰⁶Pb zircon dates (Table 2) (Tucker et al. 1990; Cohen et al. 2022).

4.c. Dicellograptus anceps Biozone

The Dicellograptus anceps biozone located in the Upper Hartfell Shale Formation constrains the end of the Late Ordovician Katian Stage (Fig. 1b) (Merriman & Roberts, 1990; Gradstein et al. 2020). A zircon fission-track date of 434 ± 12 Ma was first reported for the Dicellograptus anceps zone in Dob’s Linn (Ross, 1984). Subsequently, three (of four) multigrain zircon fractions (n = 58 of 73 grains analyzed) zircons fractions from the Dicellograptus anceps zone, located 4.5 metres below the Ordovician–Silurian boundary, yielded a mechanically (air) abraded zircon ²³⁸U-²⁰⁶Pb ID-TIMS date of 445.7 ± 2.4 Ma (Tucker et al. 1990). Using Tucker et al. (1990) zircon multigrain fractions and via spline fitting modelling, Schmitz & Ogg (2020) recalculated the age of the biozone to 444.88 ± 1.17 Ma (Schmitz & Ogg, 2020). Elsewhere, the Metabolograptus extraordinarius zone in Wanhe, SW China, the equivalent of the Dicellograptus anceps zone at Dob’s Linn, produced two CA-ID-TIMS ²³⁸U-²⁰⁶Pb dates of 443.81 ± 0.24 Ma and 444.06 ± 0.20 Ma (Table 2) (Ling et al. 2019).

5. Methods

In this study, we dated three metabentonite ash horizon samples (DL7, 19DL09 and BRS23) located within the Dicellograptus anceps, Akidograptus ascensus and Coronagraptus cyphus biozones from the Dob’s Linn biostratigraphy sections. Dob’s Linn is a Site of Special Scientific Interest (SSSI) in Scotland with restrictions on sample collection; thus, sample DL7 from the Main cliff section was provided by Richard Batchelor from archived material (Batchelor & Weir, 1988). Sample BRS23 from the Linn branch trench was provided by the British Geological Survey (Merriman & Roberts, 1990). Sample 19DL09 was collected by Catlos & Brookfield from the Linn branch trench, the same DL9 layer as in Batchelor and Weir (1988) and BRS292 in Merriman and Roberts (1990). The appropriate authorities granted permission for sample collection.

Traditional heavy mineral separation techniques were used, including deflocculation and extraction of clays via the addition of sodium hexametaphosphate and sonication to obtain maximum zircon yield. Overall, a total of 324 zircon grains were mounted in epoxy and inspected with CL using a JEOL Scanning Electron Microscope at the University of Texas at Austin, GeoMaterials Characterization and Imaging facility (GeoMatCI). Following imaging, zircons were dated using Element2 High Resolution (HR)-LA-ICP-MS in the Geo-thermochronology lab at the University of Texas at Austin. The instrument uses an Excimer (192 nm) laser ablation system and obtains isotopic measurements using ion counting. A dry ablated aerosol is introduced to the instrument by a pure He carrier gas containing the desired isotopic analytes, which for this study consist of ²³⁸U, ²³⁵U, ²³²Th, ²⁰⁶Pb, ²⁰⁷Pb and ²⁰⁸Pb. Each analysis consisted of a 2-pulse cleaning ablation, a background measurement taken with the laser off, a 30-second measurement with the laser firing and a 30 second cleaning cycle. The laser beam was 15 µm in diameter to limit analyses to specific CL domains within the zircon crystals and allow for multiple spots per grain in some cases. Elemental isotopic fractionation of Pb and Pb/U isotopes was corrected by interspersed analyses of primary and secondary zircon standards with known ages (GJ1 and Plesovice references) (Jackson et al. 2004; Sláma et al. 2008). The typical ratio of unknown standards measurements was 3:1 or 4:1. Systematic uncertainties resulting from calibration corrections are usually 1–2% for ²⁰⁶Pb/²⁰⁷Pb and ²⁰⁶Pb/²³⁸U. Pb values are reported as total Pb without any correction for potential common ²⁰⁴Pb due to isobaric interferences with ²⁰⁴Hg. Iolite software was used to process and reduce data analyses, correct instrument drift, and downhole fractionation (https://iolite-software.com/).

After LA-ICP-MS analysis, subsets of zircons from samples 19DL09 and BRS23 were removed from epoxy and subjected to CA-ID-TIMS analyses in the geochronology lab at the University of Wyoming adapted from the method of Mattinson (2005). Zircons chosen for this treatment included some of the youngest grains in BRS23 to test whether these dates reflected Pb loss and some of the oldest grains in 19DL09 to test whether these dates reflected matrix mismatch. In the CA process, zircon grains were annealed for 50 hours at 850 °C to repair fission tracks and other minor radiation damage. Zircons were then chemically abraded with HF and HNO₃ acids for 12 hours at 180 °C to partially dissolve and remove metamict portions of the grain that have experienced Pb loss due to substantial radiation damage. Single zircon grains were then spiked with a mixed ²⁰⁵Pb/²³³U/²³⁵U EARTHTIME tracer solution (ET535), dissolved in HF and HNO₃ at 235 °C for 30 hours, and converted to chlorides at 180 °C for 16 hours. Dissolved zircon samples were loaded onto single rhenium filaments with silica gel and H₃PO₄ without any further chemical processing except for three larger grains from which the Pb and U were purified on HCl-H₂O ion exchange column following Krogh (1973). Isotopic compositions were measured on a Micromass Sector 54 mass spectrometer in single collector, peak switching mode using the Daly photomultiplier collector for all isotopes (Anderson et al. 2013; Barnes et al. 2021).

Statistical values and figures, including Concordia diagrams, Kernel density estimates, Tuffzirc dates, and weighted mean distribution plots, were produced by Isoplot, Densityplotter and detritalPy (Ludwig & Mundil, 2002; Ludwig, 2008; Vermeesch, 2012; Sharman et al. 2018). A ²⁰⁶Pb/²³⁸U vs ²⁰⁷Pb/²³⁵U 10% discordance filter was implemented for all LA-ICP-MS zircon dates. Robust CA-ID-TIMS WM dates are calculated from a cluster of four or more of the youngest zircon dates overlapping within uncertainty. The youngest single grain (YSG) MDA approach is calculated from the youngest zircon date (Ludwig & Mundil, 2002). The weighted mean date (WM) is calculated from all individual zircon dates per sample using Isoplot (Ludwig, 2008). TuffZirc date is calculated using Ludwig and Mundil (2002)’s algorithm calculating the median U-Pb date of the largest coherent group of zircons dates with 2σ uncertainty using Isoplot (Ludwig, 2008). The youngest cluster of 2+ grains (YC2σ+2) is calculated from the weighted mean of the youngest zircon grain cluster of two or more grains overlapping at 2σ uncertainty (Dickinson & Gehrels, 2009). The youngest mode kernel density estimate (YMKDE) (also recognized as YPP by Dickinson and Gehrels (2009) is calculated using Vermeesch (2012)’s Densityplotter from the youngest age peak on a kernel density estimate plot (bandwidth of 10) designed from various U-Pb zircon dates while omitting single grain age peaks (Herriott et al. 2019). The youngest statistical population (YSP) is the weighted mean of the youngest subsample of two or more grains that produce a mean square weighted deviation (MSWD) close to 1 (Coutts et al. 2019). The youngest mode weighted mean (YMWM) is calculated after Tian et al. (2022), using the LA-ICP-MS zircon dates that comprise the youngest age mode from a KDE peak as a weighted mean of more than three grain overlapping at 2σ uncertainty with an approximate MSWD of 1. The KDE peak age serves as the initial reference point, with individual zircon dates extracted from both sides of the crest to attain an MSWD of 1 or an approximate value (Tian et al. 2022). The Maximum Likelihood Age (MLA) is computed via a regression algorithm employing error correlations and analytical uncertainties, assuming that data scatter primarily arises from analytical uncertainties. In the case of a correct assumption, the MSWD value should approach one (Vermeesch, 2018, 2021).

6. Results

6.a. Coronagraptus cyphus Biozone (Sample BRS23)

Samples BRS23 of the Coronagraptus cyphus zone yielded 137 zircon grains analyzed by U-Pb LA-ICP-MS and subset of 15 single grains by CA-ID-TIMS (Fig. 1b). After applying a ≤ 10% discordance filter, 133 grains ranging from Ordovician to Devonian in age were utilized to constrain an MDA for the biozone with various methods to constrain depositional ages. The youngest estimate using U-Pb LA-ICP-MS for sample BRS23 is the YSG date of 392 ± 10 Ma (5% disc), whereas the YC2σ+2 yields a date of 397 ± 10 (n = 4, MSWD = 1.40). The WM presents a date of 439 ± 2 Ma (n = 133, MSWD = 5.60), the YMKDE yields a 441 Ma date and the TuffZirc date is 441+2/−3 Ma (n = 133). The MLA produces a date of 440 ± 2 Ma (n = 133, MSWD = 5.40), and both the YSP and YMWM date is 440 ± 1 Ma (n = 83, MSWD = 1.00). CA-ID-TIMS analyses from the youngest zircon grains yielded a ²³⁸U-²⁰⁶Pb weighted mean age of 440.44 ± 0.55/0.56/0.72 Ma (±analytical/with tracer/with U-decay constant), (95% conf., MSWD 0.26, 4 of 15 analyses) (Fig. 5a; supplementary tables S1, S2).

Figure 5. Individual sample ²³⁸U-²⁰⁶Pb LA-ICP-MS dates MDA approach results: youngest single grain (YSG) (Ludwig & Mundil, 2002); youngest cluster 2+ grains at 2σ overlap (YC2 σ+2) (Dickinson & Gehrels, 2009); LA-ICP-MS total weighted mean (WM: black line); TuffZirc date (beige line) (Ludwig & Mundil, 2002); Maximum Likelihood Age (MLA) (Vermeesch, 2021); youngest mode kernel density estimate (YMKDE: pink line) (Herriott et al. 2019); youngest statistical population (YSP: yellow or brown bars when applicable) (Coutts et al. 2019); youngest mode weighted mean (YMWM: yellow bar) (Tian et al. 2022) compared to ²³⁸U-²⁰⁶Pb CA-ID-TIMS weighted mean date. (a) Sample BRS23: Coronagraptus cyphus zone. (b) Sample 19DL09: Akidograptus ascensus zone. (c) Sample DL7: Dicellograptus anceps zone.

6.b. Akidograptus ascensus Biozone (Sample 19DL09)

Sample 19DL09 of the Akidograptus ascensus zone yielded a total of 19 zircon grains analyzed by U-Pb LA-ICP-MS and subset of single grains by CA-ID-TIMS (Fig. 1b). After applying a ≤ 10% discordance filter, 17 zircon grains ranging from Ordovician to Middle Carboniferous age are utilized to constrain an MDA for the biozone with several methods to constrain depositional ages. The youngest estimate using U-Pb LA-ICP-MS for sample 19DL09 is the YSG date of 327 ± 5 Ma (1% disc), whereas the YC2σ+2 yields a date of 329 ± 13 (n = 3, MSWD = 1.70). The WM presents a date of 426 ± 22 Ma (n = 17, MSWD = 134.00), the YMKDE yields a 331 Ma date and the TuffZirc date is 447+7/−8 Ma (n = 13). The MLA produces a date of 423 ± 23 Ma (n = 17, MSWD = 110.00), and the YSP produced a date of 328 ± 5 (n = 2, MSWD = 0.92) in addition to a YMWM date of 441 ± 3 (n = 6, MSWD = 0.96). An ID-TIMS analysis without chemical abrasion from one of the youngest LA-ICP-MS dated zircon grains yielded a ²³⁸U-²⁰⁶Pb date of 339.64 ± 0.62 Ma. CA-ID-TIMS analyses produced four individual ²³⁸U-²⁰⁶Pb zircon dates from one young grain and three older plateau population grains (older recurring dates overlapping with 2σ uncertainty) previously dated by LA-ICP-MS yielding ²³⁸U-²⁰⁶Pb CA-ID-TIMS dates of 448.38 ± 1.10 Ma, 449.08 ± 1.20 Ma, 452.43 ± 3.00 Ma and 494.91 ± 1.40 Ma (Fig. 5b; supplementary tables S1, S2).

6.c. Dicellograptus anceps Biozone (Sample DL7)

Sample DL7 of the Dicellograptus anceps zone yielded a total of 40 zircon grains only analyzed by U-Pb LA-ICP-MS. After applying a ≤ 10% discordance filter, 26 grains out of 40 ranging from Ordovician to Devonian age were utilized to constrain an MDA for the biozone (Fig. 5c; supplementary table S1). The youngest estimate using U-Pb LA-ICP-MS for sample DL7 is the YSG date of 402 ± 12 Ma (5% disc), whereas the YC2σ+2 yields a date of 423 ± 9 (n = 4, MSWD = 2.00). The WM presents a date of 436 ± 4Ma (n = 26, MSWD = 7.10), the YMKDE yields a 434 Ma date and the TuffZirc date is 435+5/−2 Ma (n = 25). The MLA produces a date of 436 ± 5 Ma (n = 26, MSWD = 7.20). Both the YSP and YMWM produce a date of 433 ± 2 Ma (n = 17, MSWD = 1.00) (Fig. 5c).

7. Discussion

This study aims to re-assess the current interpretation of Dob’s Linn as the ‘GSSP’ due to the implications of understanding biological, climatic and environmental events during the Early Paleozoic. Additionally, our study is a benchmark to assess appropriate dating approaches to generate accurate MDAs of Early Paleozoic sections previously calibrated with multi-grain, ID-TIMS zircon U-Pb dates (Tucker et al. 1990; Schmitz & Ogg, 2012; Ogg et al. 2016; Gradstein et al. 2020; Cohen et al. 2022). Our study incorporates the preliminary screening of single zircon grains with CL imaging and LA-ICP-MS analyses to target the youngest and plateau populations of volcanic grains with single-grain CA-ID-TIMS analyses. This procedure permits the analysis of autocryst grain populations and filters antecrystic and/or xenocrystic zircon while mitigating the effects of Pb loss.

Concerns remain over the selection of Dob’s Linn as the global Ordovician–Silurian boundary stratotype section. According to previous studies, Dob’s Linn does not meet the international standards for a GSSP as a result of a complex tectonic and thermal history of the area affecting the stratigraphic position and accuracy of graptolite zone distributions biasing geochronology and chemostratigraphic analyses (Berry, 1987; Lesperance et al. 1987; Williams, 1988). The ICS requires a geologic section to fulfill a set of criteria to be considered a GSSP. A GSSP boundary is required to be research accessible and free to access in addition to being extensive enough to allow continuous sample collection for domestic and international researchers. A GSSP must contain a stratigraphic marker that defines the lower boundary of a geologic Stage. The boundary must present diversity and abundance of well-preserved fossils throughout the boundary interval, including secondary markers such as other fossils and chemical changes manifested in regional and global stratigraphic sections. The stratigraphic section must have layers containing minerals that can be radiometrically dated and adequate thickness allowing global correlation, including continuous sedimentation without gaps or changes in facies. The boundary is required to be unaffected by tectonic disturbances and metamorphism (Remane et al. 1996; Gradstein & Ogg, 2020).

Initially, Dob’s Linn was selected in 1979 by the Boundary Working Group as the GSSP for the base of the Parakidograptus acuminatus zone marking the base of the Silurian and thus the Ordovician–Silurian boundary and later reassessed to the Akidograptus ascensus graptolite zone (Fig. 1b) (Ross, 1984; Cocks, 1985, 1988; Rong et al. 2008). The primary concerns for questioning Dob’s Linn as a reference section is due to the limited lateral extent of graptolite zones, scarcity of fossils other than graptolites, and the locality’s tectonic and thermal disturbed sections forming large and micro-scale folds and faults across the Moffat Shale, disputing the accuracy of the graptolite data (Leggett et al. 1979; Williams, 1983; Williams & Rickards, 1984; Berry, 1987; Lesperance et al. 1987; Williams, 1988). Isotopic carbon data points to the Metabolograptus persculptus zone as a possibility this graptolite horizon can be used as the Ordovician–Silurian boundary rather than the current assessed Akidograptus ascensus zone (Fig. 1b) (Berry, 1987). Additionally, the Ordovician–Silurian stratigraphic section from Anhui, China, is reported to have an ideal abundance and diversity of graptolites without tectonic disturbances making it an ideal candidate for the Ordovician–Silurian boundary GSSP (Ji-jin et al. 1984; Berry, 1987).

Due to Dob’s Linn’s SSSI status, it proved difficult to collect and obtain adequate sample sizes to generate the large quantities of zircon ideally used to produce robust ²³⁸U-²⁰⁶Pb dates (Vermeesch, 2004; Andersen, 2005). However, with the samples provided, we were able to generate meaningful results. In addition, a comparison between LA-ICP-MS and CA-ID-TIMS results for the same grains provides some important observations that should be made when assessing MDAs using the laser-based approach alone.

In the case of this study, the youngest chronostratigraphic sample is BRS23 from the Coronagraptus cyphus zone, presenting an LA-ICP-MS KDE distribution with a primary peak of 441 Ma showing a younger skewed tail incorporating Devonian zircon dates as young as 392 ± 10 Ma to as old as Ordovician 484 ± 13 Ma (Fig. 5a). The range of LA-ICP-MS dates from sample BRS23 are significantly younger and older than its currently recognized Silurian age of 439.57 ± 1.33 Ma (Tucker et al. 1990; Gradstein et al. 2020; Schmitz & Ogg, 2020), and this study’s CA-ID-TIMS ²³⁸U-²⁰⁶Pb WM date of 440.44 ± 0.72 Ma The older CA-ID-TIMS dates from sample BRS23 confirm the presence of antecrysts with pre-eruptive growth in the zircon grains within this metabentonite (Wotzlaw et al. 2013; Schaltegger et al. 2014). As shown in Fig. 6, zircon grains with muted zoning textures have dates that are indicative of ash fall origins (autocrysts), and grains with oscillatory zoning as a result of episodic magmatic growths tend to be associated with antecrysts. The muted CL may reflect rapid crystallization of eruptive zircons from a single homogenous magma and could be useful in differentiating them from antecrysts. Due to detecting both autocrysts and antecrysts in this single bentonite layer at Dob’s Linn, caution is necessary when dating this section with ID-TIMS multigrain zircon fractions without any zircon grain pre-screening by CL imaging or LA-ICP-MS.

Figure 6. CL images of representative Dob’s Linn metabentonite zircon grains with ²³⁸U-²⁰⁶Pb LA-ICP-MS and CA-ID-TIMS dates. Ovals indicate locations of LA-ICP-MS analyses. Zircons with broad, muted zoning textures may reflect eruptive zircons, whereas tighter, concentric, oscillatory zoning is more typical of magmatic zircon growths. Bright cathodoluminescence zones correspond to high U content.

Furthermore, it is important to note that most of the single-grain comparative dates between LA-ICP-MS and CA-ID-TIMS overlap within uncertainty except for the youngest and in some cases the oldest LA-ICP-MS dates (Fig. 7). For sample BRS23, Figs. 6 and 7 compares LA-ICP-MS and CA-ID-TIMS dates from the same individual grains where two of the youngest LA-ICP-MS population grains show differences of up to 40 Ma with CA-ID-TIMS showing Pb loss is present and effectively removed by the chemical abrasion treatment. In addition, LA-ICP-MS dates that are older than their CA-ID-TIMS dates (Figs 6 and 7; supplementary table S2) may reflect a mismatch in ablation rates between samples and standards that lead to bias in the U-Pb downhole fractionations. We refer to this matrix mismatch as it likely stems from different crystal lattice states of samples and standards. Although this effect can theoretically produce dates that are both too young and too old, standards are typically low in U and have less lattice damage than many samples, so the effect is skewed towards under-representation of U or apparent U loss and dates that are too old. The Concordia diagram in Fig. 8a used to evaluate the age consistency between the two chronometers ²³⁸U-²⁰⁶Pb and ²³⁵U-²⁰⁷Pb and disturbances within the U-Pb system by Pb loss displays a prominent age cluster. However, the BRS23 Concordia diagram also shows younger than expected clusters of LA-ICP-MS dates exhibiting Pb loss in the system. Using various LA-ICP-MS MDA calculation methods for sample BRS23, the YSG date of 392 ± 10 Ma (5% disc) and the YC2σ+2 with a date of 397 ± 10 (n = 4, MSWD = 1.40) produced the youngest MDA dates for this sample (Fig. 5a). The WM, MLA and TuffZirc date yielded dates of 439 ± 2 Ma (n = 133, MSWD = 5.60), 440 ± 2 Ma (N = 133, MSWD = 5.40) and 441+2/−3 Ma (n = 133) overlapping within a larger uncertainty with the current interpreted age of 439.57 ± 1.33 Ma by Schmitz & Ogg (2020), and this study’s CA-ID-TIMS ²³⁸U-²⁰⁶Pb weighted mean date of 440.44±0.55/0.56/0.72 Ma (±analytical/with tracer/with U-decay constant) (95% conf., MSWD 0.26, 4 of 15 analyses). However, the YMKDE with a date of 441 Ma and both the YSP and YMWM with the same date of 440 ± 1 Ma (n = 83, MSWD = 1.00) approximate the current interpreted age of the biozone and our CA-ID-TIMS date with higher precision compared to the WM and Tuffzirc dates (Figs. 5a, and 9).

Figure 7. One-to-one comparison between individual LA-ICP-MS U-Pb (blue) vs. CA-ID-TIMS (red) U-Pb dates. Error bars represent 2σ uncertainty.

Figure 8. U-Pb LA-ICP-MS and CA-ID-TIMS Concordia diagram reported as total Pb for zircon data from metabentonites in the Hartfell and Birkhill shales at Dob’s Linn, Scotland. (a) Sample BRS23: Coronagraptus cyphus zone, accepted age from Tucker et al. (1990). (b) Sample 19DL09: Akidograptus ascensus zone, accepted age from Cohen et al. (2022). (c) Sample DL7: Dicellograptus anceps zone, accepted age from Tucker et al. (1990). Black ovals show LA-ICP-MS dates, and red ovals display CA-ID-TIMS dates.

Figure 9. Hartfell and Birkhill Shale stratigraphy with graptolite biozones from Dob’s Linn (after Batchelor & Weir, 1988; Merriman & Roberts, 1990). Samples BRS23 and 19DL09 from the Linn trench branch and DL7 from the Main cliff location. Comparison between currently recognized zone ages and new U-Pb dates presented in this study. D. anceps zone age from Tucker et al. (1990); A. ascensus zone age defining Ordovician–Silurian boundary (dashed line) from Cohen et al. (2022); C. cyphus zone age from Tucker et al. (1990). LA-ICP-MS YMWM and TuffZirc date MDA approaches present suitable dates for samples BRS23 and 19DL09 to currently recognized biozone ages, and this study’s BRS23 CA-ID-TIMS date. Sample DL7 displays potential stratigraphic misplacement or significant Pb loss presented by the younger MDA dates than sample BRS23.

The Akidograptus ascensus zone associated with sample 19DL09 in this study is interpreted as the global standard reference section for the Ordovician–Silurian boundary with a calculated age of 443.8 ± 1.5 Ma with the use of spline fitting interpolation using various radioisotopic dates (Rong et al. 2008; Cohen et al. 2022). Sample 19DL09 in this study presents the first ²³⁸U-²⁰⁶Pb zircon dates from a Dob’s Linn metabentonite in the Akidograptus ascensus zone. The LA-ICP-MS KDE distribution shows two bimodal peaks, including Carboniferous zircon dates as young as 327 ± 5 Ma to as old as Ordovician 464 ± 7 Ma (Fig. 5b). Three significantly young Carboniferous zircon dates form the youngest peak shown in the KDE distribution. The second older peak comprises 14 zircon dates, from which most are Ordovician–Silurian age, except for one young Silurian–Devonian date. The range of LA-ICP-MS dates from sample 19DL09 is predominantly skewed towards significantly younger zircon dates than the current interpreted Ordovician–Silurian boundary age of 443.8 ± 1.5 Ma (Cohen et al. 2022). Using various LA-ICP-MS MDA calculation methods for sample 19DL09, the youngest MDA dates were produced with the YSG, YC2σ+2, WM, YMKDE and YSP approaches (Fig. 5b). The YSG yielded a date of 327 ± 5 Ma (1% disc), and the YC2σ+2 produced a date of 329 ± 13 (n = 3, MSWD = 1.70). The WM yielded a date of 426 ± 22 Ma (n = 17, MSWD = 134.00) and the MLA produced a date of 423 ± 23 Ma (n = 17, MSWD = 110.00). The YMKDE shows a date of 331 Ma, and the YSP produced a date of 328 ± 5 (n = 2, MSWD = 0.92). Only the TuffZirc date with a date of 447+7/−8 Ma (n = 13) and the YMWM with a date of 441 ± 3 (n = 6, MSWD = 0.96) calculated the current Ordovician–Silurian boundary age of 443.8 ± 1.5 Ma within uncertainty (Fig. 9).

Due to the destructive nature of the CA-ID-TIMS method, not all the youngest LA-ICP-MS dated zircon grains endure the chemical abrasion process, thus preventing this study from producing a robust ²³⁸U-²⁰⁶Pb CA-ID-TIMS WM date for sample 19DL09. However, we present four individual CA-ID-TIMS zircon dates previously screened with LA-ICP-MS from a young and three older plateau population grains with ²³⁸U-²⁰⁶Pb CA-ID-TIMS dates of 448.38 ± 1.10 Ma, 449.08 ± 1.20 Ma, 452.43 ± 3.00 Ma and 494.91 ± 1.40 Ma (supplementary table S2). As shown in Fig. 5, the zircon grains from sample 19DL09 also display muted zoning textures reflecting ash fall origins (autocrysts) and oscillatory textured grain from previous magmatic growths (antecrysts) (Wotzlaw et al. 2013; Schaltegger et al. 2014). Two of the CA-ID-TIMS analyses, when compared to their respective LA-ICP-MS dates, yield slightly younger zircon dates with overlapping uncertainty (Fig. 7b). One grain shows a CA-ID-TIMS date slightly younger by 2 Ma than its LA-ICP-MS date (Figs. 6, and 7b). However, one anomalous grain from the youngest population of LA-ICP-MS dates produced a significantly older CA-ID-TIMS date by 160 Ma, yielding a CA-ID-TIMS date older than the stratigraphic age representing the inclusion of an inherited core based on its discordance (sample 19DL09 g2; Figs. 7b and 8b; supplementary table S2). The differences between the LA-ICP-MS and CA-ID-TIMS individual grain comparison for sample 19DL09 are interpreted to reflect both matrix mismatches and Pb loss effects for these zircon grains due to the younger and older dates when compared to CA-ID-TIMS dates. Concordia diagram in Fig. 8b illustrates two distinguishable clusters with the youngest cluster of LA-ICP-MS date manifesting significant Pb loss. The significant Pb loss effect is also reflected with the younger LA-ICP-MS younger zircon population (Fig. 5b). Furthermore, one grain from the youngest population of LA-ICP-MS dates with an anomalous Carboniferous young date of 338 ± 12 was not chemically abraded and dated with ID-TIMS yielding the same unusually young date of 339.64 ± 0.62 (supplementary table S2). These zircon dates yielding the same young dates support our interpretation that Pb loss is the dominant factor with the anomalously young grains dated using LA-ICP-MS and the chemical abrasion process is effective at removing the effects of Pb loss.

The oldest chronostratigraphic sample is DL7 from the Dicellograptus anceps zone showing a symmetric LA-ICP-MS KDE distribution with a single broad peak of 434 Ma incorporating Devonian zircon dates as young as 402 ± 12 Ma to as old as Ordovician 462 ± 10 Ma (Fig. 5c). The range of LA-ICP-MS dates from sample DL7 visually does not look skewed towards younger dates, but the majority of the grains yield predominantly younger zircon dates than its biozone’s interpreted age of 444.88 ± 1.17 Ma (Gradstein et al. 2020; Schmitz & Ogg, 2020). Although we were unable to produce CA-ID-TIMS dates for sample DL7 due to loss of zircons to complete dissolution during the chemical abrasion process, we present several LA-ICP-MS MDA calculation methods for the Dicellograptus anceps zone at Dob’s Linn. The abundance of younger grains along with the complete dissolution of them is consistent with metamict zircons. All MDA approaches for sample DL7, including the YSG, YC2σ+2, WM, TuffZirc date, MLA, YMKDE, YSP and YMWM, yielded significantly younger or inconsistent dates than the current assessed biozone age of 444.88 ± 1.17 Ma (Fig. 5c). The YSG yielded a date of 402 ± 12 Ma (5% disc), and the YC2σ+2 produced a date of 423 ± 9 (n = 4, MSWD = 2.00). The WM yielded a date of 436 ± 4 Ma (n = 26, MSWD = 7.10), and the TuffZirc date produced a date of 435+5/−2Ma (n = 25). The MLA produces a date of 436 ± 5 Ma (n = 26, MSWD = 7.20). The YMKDE yielded a date of 434 Ma, and both the YSP and YMWM yielded the same date of 433 ± 2 (n = 17, MSWD = 1.00) (Figs. 5c, and 9). It is important to note that sample DL7, linked to the Dicellograptus anceps zone at Dob’s Linn, comes from the Main Cliff locality rather than the Linn Branch GSSP location, only separated by several hundred metres horizontally (Batchelor & Weir, 1988; Williams, 1988; Verniers & Vandenbroucke, 2006). The Concordia diagram in Fig. 8c shows dispersed zircon dates including a minor cluster with large uncertainties approximating the current interpret age of 444.88 ± 1.17 Ma. However, a multitude of dates are significantly younger than expected manifesting significant Pb loss. The significant Pb loss effect is also reflected in the overall LA-ICP-MS KDE zircon distribution (Fig. 5c). Furthermore, one grain from the youngest population of LA-ICP-MS dates with an anomalous Carboniferous young date of 338 ± 12 was not chemically abraded and dated with ID-TIMS yielding the same unusual young date of 339.64 ± 0.62 (supplementary table S2). Based on the abundance of young LA-ICP-MS dates from sample DL7 showing younger MDA calculations than the youngest BRS23 sample in this study, the zircon grains in this particular horizon may have experienced more significant Pb loss, or the Dicellograptus anceps zone is incorrectly assigned in the Main Cliff locality (Fig. 9). As previously mentioned, the Dob’s Linn locality is considerably tectonically and thermally disturbed where the same graptolite biozone identification occurs in different stratigraphic positions separated by large- and small-scale faults (Berry, 1987; Lesperance et al. 1987). Graptolite horizons in the Main Cliff locality can potentially be misplaced in the stratigraphy, thus not presenting the first appearance of a specific graptolite fauna but a later occurrence.

All three samples in this study demonstrate that MDA interpretations using YSG and YC2σ +2 are invalid for this Early Paleozoic section. These approaches can be considered less conservative and follow the theory of using the youngest concordant zircons or youngest concordant zircon clusters as the maximum age of an enclosing sediment, thus often generating considerably younger dates than the true depositional age (TDA) (Herriott et al. 2019). As shown with this study’s YSG and YC2σ +2 dates, discordance alone is not adequate to identify Pb loss from Phanerozoic LA-ICP-MS data (Anderson et al. 2019). In all three samples (BRS23, 19DL09, DL7), the LA-ICP-MS YSG and YC2σ +2 dates produced anomalous young dates (Fig. 5). For sample DL7, there are two possibilities (or possibly a combination) as to why we identify complications with the sample. The zircons experienced more Pb loss generating biases with the calculated MDA approaches or incorrect stratigraphic assignment. In the case of sample BRS23, the YSG and YC2σ +2 yield much younger dates by 36 Ma and 31 Ma when compared with our CA-ID-TIMS date of 440.44 ± 0.72 Ma and the current recognized Coronagraptus cyphus zone age of 439.57 ± 1.33 Ma (Gradstein et al. 2020; Schmitz & Ogg, 2020). For sample 19DL09, the YSG and YC2σ +2 yield significantly younger dates by 110 Ma and 100 Ma compared to the current recognized Akidograptus ascensus age of 443.8 ± 1.5 Ma (Cohen et al. 2022). For sample DL7, the YSG and YC2σ +2 yield younger dates by 30 Ma and 12 Ma compared to the current recognized Dicellograptus anceps age of 444.88 ± 1.17 Ma (Gradstein et al. 2020; Schmitz & Ogg, 2020). The YMKDE yielded a suitable date of 441 Ma for sample BRS23 within uncertainty of its current assessed age. However, the YMKDE produced considerably younger dates of 110 Ma and 10 Ma for samples 19DL09 and DL7. The YSP approach only produced a suitable date for sample BRS23, comparable to this study’s CA-ID-TIMS date and current recognized biozone age within uncertainty. However, the YSP yielded younger dates by 112 Ma and 9 Ma for samples 19DL09 and DL7. The TuffZirc date yielded appropriate dates for samples BRS23 and 19DL09, though with a considerably larger uncertainty when compared to our CA-ID-TIMS date or current recognized biozone age. In the case of sample DL7, the TuffZirc date provided a younger date by 4 Ma. Similarly, the YMWM also produced suitable dates for samples BRS23 and 19DL09 with improved uncertainty compared to the TuffZirc date, and within uncertainty of our CA-ID-TIMS date or current recognized biozone age. However, for sample DL7, the YMWM yields a date 9 Ma younger than its current assessed age.

The study’s results suggest that both LA-ICP-MS and CA-ID-TIMS dating approaches are needed in localities like Dob’s Linn, where extensive post-depositional structural and hydrothermal alterations produce a high percentage of discordant and potentially metamict zircons due to the widespread tectonic activity associated with the formation of the Caledonian mountains (Fig. 4) (Lesperance et al. 1987; Chew & Strachan, 2014). Integrating LA-ICP-MS and CA-ID-TIMS provides the benefit of pre-screening and eliminating possible older detrital grains and identifying target zircons from the youngest populations for CA-ID-TIMS analyses. The complex tectonic activity at Dob’s Linn distorted the graptolite biozones in southern Scotland, inducing biostratigraphic misrepresentation and potentially influencing an increase in metamict grains, thus overall inducing younger U-Pb zircon dates due to Pb loss (Lesperance et al. 1987).

Comparative single grain dates between LA-ICP-MS and CA-ID-TIMS overlap within uncertainty predominantly with zircon plateau populations. However, this is not the case for the youngest individual grains and clusters with differences of more than 100 Ma due to Pb loss and matrix mismatch influence. Additionally, individual LA-ICP-MS dates that are older than CA-ID-TIMS dates are interpreted to reflect a matrix mismatch that biases the U-Pb downhole fractionation (Figs. 5, and 7, supplementary table S2). The best LA-ICP-MS MDA estimates are generated from calculations using averages and not primarily the youngest grains so that both Pb loss and matrix mismatch effects are minimized (e.g., Tian et al. 2022). In cases where significant number of zircons have experienced Pb loss, the chemical abrasion process becomes crucial, as averages may not provide accurate results. Our study’s findings support the conclusions of Tian et al. (2022) with the TuffZirc date and YMWM yielded the best results to achieve appropriate MDAs estimation when only LA-ICP-MS data is available. For this study, the LA-ICP-MS MDA methodologies such as the TuffZirc date and YMWM yielded results that were in line with CA-ID-TIMS analysis, within the margin of uncertainty. This outcome was achieved by employing average calculations that encompassed both the youngest grains affected by variable Pb loss and moderately older grains that remained unaffected by significant Pb loss. The TuffZirc date by Ludwig and Mundil (2002) demonstrated suitable results, although with larger uncertainty for two of our three samples, and generated a young MDA for the third complex sample (DL7) by only 4 Ma younger than its current assessed age (Fig. 9). The YMWM from Tian et al. (2022) also generated appropriate results, though more robust with superior uncertainty than the TuffZirc date for two of our three samples; however, for the third problematic sample (DL7), the approach determined an MDA date 9 Ma younger than its current assessed age (Fig. 9). The TuffZirc date and YMWM present younger dates for the older Dicellograptus anceps zone (sample DL7) than the current recognized age and our CA-ID-TIMS date for the Coronagraptus cyphus zone (sample BRS23), which is considered the youngest stratigraphic sample in this study.

Based on the GSSP requirements by the ICS, including the findings from previous studies describing the inconsistencies of Dob’s Linn as a reference section, and our finding in this study, such as the potential of the Dicellograptus anceps zone being incorrectly assigned stratigraphically, Dob’s Linn’s GSSP status appears to be questionable (Berry, 1987; Lesperance et al. 1987; Remane et al. 1996). We suggest the re-examination of Dob’s Linn, both with biostratigraphy and additional sample collection for future CA-ID-TIMS zircon dates to improve accuracy and precision or considering other Ordovician–Silurian boundary outcrops such as the ones in Anticosti Island, Canada or South China for future studies involving the Ordovician–Silurian periods.

8. Conclusions

This study presents new data from the Dicellograptus anceps, Akidograptus ascensus and Coronagraptus cyphus zones at Dob’s Linn, Scotland. We produced a high-precision CA-ID-TIMS date of 440.44 ± 0.55/0.56/0.72 Ma (±analytical/with tracer/with U-decay constant) for the Coronagraptus cyphus zone. Comparisons between CA-ID-TIMS and LA-ICP-MS U-Pb zircon dates for the metabentonites encompassing the Akidograptus ascensus and Coronagraptus cyphus zones demonstrate the presence of both autocrysts and antecrysts in addition to significant Pb loss and matrix mismatch between LA-ICP-MS unknowns and standards (Figs. 5–7; supplementary table S2). The presence of autocrysts and antecrysts in Dob’s Linn holds significance due to its impact on the accuracy of the established biozone ages. Previously, these ages were determined with multigrain zircon fractions ID-TIMS analyses (Tucker et al. 1990). Comparative single U-Pb dates between LA-ICP-MS and CA-ID-TIMS overlap within uncertainty primarily with zircon plateau populations; however, this is not the case for the youngest grains and youngest cluster populations showing anomalous differences of more than 100 Ma with the currently assessed biozone ages and our CA-ID-TIMS dates (supplementary table S2). The Pb loss and matrix mismatch is corroborated with the notably younger zircon dates and older individual LA-ICP-MS dates compared to individual CA-ID-TIMS analyses (Figs. 6 and 7; supplementary table S2).

We suggest integrating LA-ICP-MS and CA-ID-TIMS whenever possible for MDA calculations to screen and eliminate older detrital grains and focus on the youngest individual grains and populations for CA-ID-TIMS analyses. In cases where CA-ID-TIMS analysis is not feasible, we strongly advocate the annealing of both unknown and standard zircon gains to enhance and standardize matrix conditions (Allen & Campbell, 2012). MDAs based on a small number of grains (i.e., YSG, YC2σ +2) are unreliable in our study. We recommend utilizing MDA calculations by considering averages of grains beyond solely relying on the youngest zircon grains to mitigate potential issues related to Pb loss and matrix mismatch effects. MDA methodologies such as the TuffZirc date and YMWM demonstrated optimal performance attributed to the incorporation of older co-genetic LA-ICP-MS zircon dates that remained unaffected by substantial Pb loss. As a result, these older co-genetic dates superseded the influence of younger grains impacted by variations of Pb loss. Based on our results, the TuffZirc date and YMWM produced adequate MDA calculations when only LA-ICP-MS data is available as they yield comparable results to our CA-ID-TIMS analyses or the currently recognized biozone ages within uncertainty (Fig. 9) (Ludwig & Mundil, 2002; Tian et al. 2022). In the case of sample DL7, we are uncertain whether sample DL7 associated with the Dicellograptus anceps zone at Dob’s Linn reflect Pb loss, stratigraphic misplacement, or both due to widespread tectonic and thermal activity. More sample material from the Dob’s Linn locality is necessary to acquire additional CA-ID-TIMS analyses. The LA-ICP-MS TuffZirc date and YMWM MDA approaches indicate younger dates for the Dicellograptus anceps zone than the youngest sample of the study BRS23 from the Coronagraptus cyphus zone (Fig. 9). The potential biostratigraphy and stratigraphic misplacement encountered with this study, along with the International Commission of Stratigraphy (ICS) GSSP requirements and previous reports of the inadequacy of Dob’s Linn as a global reference section, raises concerns on the validity of Dob’s Linn as the Ordovician–Silurian GSSP type section (Berry, 1987; Lesperance et al. 1987; Remane et al. 1996). A comprehensive future re-examination of Dob’s Linn is essential using biostratigraphy and geochronology to assess the legitimacy of Dob’s Linn as a GSSP or the appointment of a new proper location as the Ordovician–Silurian boundary GSSP.