Kikiktat volcanics of Arctic Alaska—Melting of harzburgitic mantle associated with the Franklin large igneous province

The record of mafic magmatism within Laurentia during the mid-Neoproterozoic (850–720 Ma) is dominated by magmatic events associated with the polyphase breakup of Rodinia (Ernst et al., 2008; Li et al., 2008). The most significant magmatic episodes include the emplacement of the ca. 825 Ma Gairdner large igneous province in Australia (Sun and Sheraton, 1996; Wingate et al., 1998; Zhao et al., 1994) and the synchronous Guibei large igneous province in China (Wang et al., 2008, 2010), the ca. 780 Ma Gunbarrel large igneous province in North America (Harlan et al., 2003) and the synchronous Kanding large igneous province in China (Zhang et al., 2009), the ca. 755 Ma Mundine Wells event of Australia (Li et al., 2006; Wingate and Giddings, 2000), and the ca. 720 Ma Franklin large igneous province in North America and Greenland (Denyszyn et al., 2009; Macdonald et al., 2010). The temporal and spatial record of this basaltic volcanism is inconsistent with a single plume-generated model for these continental flood basalts; rather, it broadly appears to record diachronous plume magmatism and supercontinent breakup (Ernst et al., 2008; Li et al., 2008).

Igneous rocks associated with these events are mostly dominated by basalt to basaltic andesites (hypersthene to quartz-normative tholeiites; Fig. 1), typical of continental flood basalts, many of which are too evolved to have been in equilibrium with peridotitic mantle (i.e., in equilibrium with olivines of less than Fo₈₉ composition). In fact, with average MgO compositions of ∼6% and SiO₂ compositions of ∼52% (Fig. 1A), the average composition is close to being a basaltic andesite; consequently, most of these continental flood basalts have undergone significant gabbroic fractionation (i.e., olivine + plagioclase + clinopyroxene; Fig. 1B), making it difficult to reconstruct a composition that would have been in equilibrium with mantle peridotite.

Another complicating factor in understanding the petrogenesis of these volcanics is that in many cases, only the feeder dikes to the now-eroded volcanics have been preserved. In these cases, it may be hard to distinguish true liquid compositions from mixtures that were modified by cumulate phases. Therefore, it is difficult to ascertain the relative contribution of the subcontinental lithospheric mantle, asthenospheric mantle, or plumes in generating these significant magmatic events.

The ca. 720 Ma Kikiktat volcanics (previously termed the Mount Copleston volcanics by Moore, 1987) of the North Slope of Arctic Alaska represent a suite of exceptionally well-preserved lavas that lie directly below 717–662 Ma Sturtian glacial deposits (Macdonald et al., 2009, 2010; Rooney et al., 2014; Strauss et al., 2013). Occurring as both subaerial flows with preserved pahoehoe flow tops as well as subaqueous pillow basalts, these lavas were first described by Moore (1987) as continental tholeiites. However, detailed petrographic and geochemical analyses presented here reveal that these volcanics are picritic melts, including high-MgO picrites that have only undergone olivine fractionation, whereas more-evolved samples are picro-basalts that lie along a gabbroic fractionation trend presenting a coherent liquid line of descent from a high-MgO (∼12%) parental melt. In this paper, we present new geochronological constraints that link these basalts to the ca. 720 Ma Franklin large igneous province. We also present new major, trace, and isotopic data for this unique suite of continental flood basalts and address the tectonic and thermal regime necessary for their petrogenesis, helping to elucidate the magmatic history of the Franklin large igneous province and the tectonic history of Arctic Alaska.

Pre-Mississippian rocks of the Shublik and Sadlerochit Mountains in the northeastern Brooks Range of Arctic Alaska are exposed in a series of north-vergent thrust sheets of the Jurassic–Paleogene Brookian fold-and-thrust belt (Fig. 2; Wallace and Hanks, 1990). These deposits form part of the North Slope subterrane of the Arctic Alaska–Chukotka terrane (cf. Strauss et al., 2013), a Cordilleran “suspect” terrane that has an ambiguous origin and displacement history (e.g., Amato et al., 2009; Miller et al., 2010, 2011). Although previous studies have suggested an exotic origin for the entire Arctic Alaska–Chukotka terrane (Amato et al., 2009; Blodgett et al., 2002; Macdonald et al., 2009), recent work hints at a more complex history that most likely involved the juxtaposition of multiple allochthonous and parautochthonous Precambrian–Paleozoic crustal fragments to form the composite Mesozoic Arctic Alaska–Chukotka terrane (Strauss et al., 2013). Since this work concentrates on the Arctic Alaska portion of the composite Arctic Alaska–Chukotka terrane, we will use Arctic Alaska in the rest of this manuscript.

The Kikiktat volcanics have been mentioned in a number of geological mapping and stratigraphic studies, referred to informally as the Mount Copleston volcanics (Blodgett et al., 1992, 2002; Clough and Goldhammer, 2000; Dutro, 1970; Macdonald et al., 2009; Moore, 1987; Robinson et al., 1989; Strauss et al., 2013). However, this name was also used to formalize an Early Devonian carbonate unit in the Shublik and Sadlerochit Mountains—the Mount Copleston Limestone, which unconformably overlies the Neoproterozoic–Ordovician Katakturuk Dolomite and Nanook Limestone (Blodgett et al., 1992). Therefore, to avoid repetition of the name, we refer to these volcanics rocks as the Kikiktat volcanics based on their spectacular exposure near Kikiktat Mountain on the Hula Hula River (Macdonald et al., 2009; Reed, 1968; Reiser, 1971).

Clough et al. (1990) determined a Rb-Sr whole-rock isochron age of 801 ± 20 Ma along with a Sm-Nd isochron age of 704 ± 38 Ma on a feeder sill to the Kikiktat volcanics from the eastern Sadlerochit Mountains. Macdonald et al. (2009) and Strauss et al. (2013) presented detrital zircon data from sandstone units underlying and interbedded with the volcanics that contain grains as young as 760 ± 11 Ma. Directly overlying these basalts, there is a discontinuous volcaniclastic conglomerate and a carbonate matrix diamictite, informally termed the Hula Hula diamictite, which Macdonald et al. (2009) suggested has lithostratigraphic and chemostratigraphic ties to widespread 717–662 Ma Sturtian-age glacial diamictites (Macdonald et al., 2009, 2011, 2012). This correlation is confirmed herein by new U-Pb zircon ages that link these volcanics temporally with the ca. 720 Ma Franklin large igneous province (Heaman et al., 1992; Macdonald et al., 2010).

The Kikiktat volcanics are variably exposed at the base of major north-vergent thrust sheets in the Shublik and Sadlerochit Mountains and in the core of a large Brookian anticline near Kikiktat Mountain (Fig. 2). The general thickness patterns and flow geometries of the Kikiktat volcanics are poorly constrained due to their common structural imbrication and limited exposure; however, one can construct a decent framework for their field relations and characteristics by incorporating observations from different exposures.

In the western Shublik Mountains, the rusty-weathering, dark-green Kikiktat volcanics are best exposed in outcrops along Nanook Creek and are at least 450 m thick (Fig. 3). Here, the basal portion of the volcanics consists of ∼80 m of structurally imbricated, heavily altered, and weakly vesiculated pillowed lavas interbedded at the top with ∼30 m of maroon and green phyllite and minor fine- to medium-grained quartz and lithic arenite that hosts ripple cross-stratification, heavy mineral lamination, and meter-scale upright isoclinal folds. Despite the structural complexities at the base of this section, we agree with Moore (1987) that the small isoclinal folds are synsedimentary and associated with the eruption of lava onto wet, partially lithified sediment, given their intimate relationship with bounding flows. This is also supported by evidence for soft-sediment deformation and flame structures in the interbedded siliciclastic strata. The contact between this siliciclastic unit and the thicker upper portion of the volcanics is recessive and not exposed (Fig. 3).

The more-continuous upper section of the Kikiktat volcanics consists of ∼350 m of stacked 0.7–32-m-thick flows with abundant pillows and evidence for subaerial exposure in the form of distinct scoriaceous breccias, well-preserved pahoehoe textures, and occasional intervals of reworked volcaniclastic rocks (Fig. 4A). Flow tops are easily recognizable by horizons of heavy vesiculation and red-weathering zones of concentrated volcanic breccia, both of which are generally capped by sharply defined bases of overlying flows. Locally, the more massive flows display distinct columnar jointing and concentrated large amygdules (0.1–0.5 cm) at flow tops that are filled with quartz, calcite, chlorite, and occasional coatings of native copper and iron oxides. Disseminated native copper is common throughout the Kikiktat volcanics (Moore, 1987).

Near Kikiktat Mountain (Fig. 4B), the volcanics are greater than 500 m thick and consist of heavily altered and chloritized metabasalts with occasional evidence for pillows and abundant intervals of volcaniclastic breccia (Macdonald et al., 2009; Reed, 1968; Reiser, 1971). In contrast, in the eastern Sadlerochit Mountains along Nularvik Creek (Fig. 4C), the Kikiktat volcanics are well preserved and host similar characteristics to the Shublik Mountains; however, they are only ∼105 m thick (Macdonald et al., 2009). It is not clear if this thinning is due to stratigraphic or structural truncation. The contact between the Kikiktat volcanics and the underlying siliciclastic strata of the Neruokpuk Formation (sensu Robinson et al., 1989), or O.G. map unit (Macdonald et al., 2009), is generally structural, although at some localities in the eastern Sadlerochit Mountains, it is most likely conformable (Macdonald et al., 2009; Robinson et al., 1989). The upper contact of the Kikiktat volcanics and the Hula Hula diamictite and/or Katakturuk Dolomite is clearly disconformable and locally erosional. For example, in the Shublik Mountains, the Hula Hula diamictite is missing, and a heavily oxidized and sharp surface marks the contact between the Katakturuk Dolomite and underlying flows of the Kikiktat volcanics (Fig. 4C). In contrast, near Kikiktat Mountain, over 50 m of volcaniclastic phyllitic schist and diamictite (Fig. 4B) rest between pillowed lavas of the Kikiktat volcanics and deep-water debris-flow deposits of the Katakturuk Dolomite. Zircons extracted from this local volcaniclastic interval (Fig. 4A) yield the new geochronological constraints presented here. Importantly, these samples are from a green-weathering volcaniclastic interval between maroon volcaniclastic rocks without evidence for glaciation and carbonate matrix diamictite with exotic clasts and matrix-supported lonestones that suggest a glacigenic origin (Macdonald et al., 2009); therefore, the geochronological data provide a synchronous to minimum age constraint on volcanism in the Kikiktat volcanics and a maximum age constraint on the Hula Hula diamictite.

The maximum metamorphic grade reached by the Kikiktat volcanics is greenschist facies, with chlorite alteration a ubiquitous feature of all samples. Serpentinization is locally present in high-MgO basalt samples, and any original groundmass is now composed of a fine-grained sericitic mineral assemblage (Fig. 5A). There is evidence for secondary silicification in some samples, and the average loss on ignition value for the sample suite is ∼3%.

High-Mg basalt samples show abundant small phenocrysts of olivine with extensive iddingsite alteration dominated by hematite (Figs. 5A–5F). Oxide dominance over hydrous iron phases and Mg clays is characteristic of the high-temperature variety of iddingsite alteration (Baker and Haggerty, 1967; Goff, 1996), providing evidence that this alteration was contemporary with eruption. While many of the olivine crystals are solid and polyhedral, many show skeletal features (Figs. 5B and 5C). Abundant plagioclase laths and subophitic interstitial clinopyroxene are the other dominant minerals in Kikiktat samples. High-Mg samples also contain spherulitic and variolitic textured clinopyroxene, the growth of which has displaced olivine grains to the outside of the spherules and provides evidence that this is a feature of crystallization and not a devitrifying feature from glass (Figs. 5D and 5E). Spherulitic/variolitic textured clinopyroxene and skeletal olivine suggest that these high-Mg lavas were superheated with respect to olivine and clinopyroxene (Fowler et al., 2002).

Low-Mg samples are dominated by plagioclase, clinopyroxene, and a glassy altered groundmass. Very rare large phenocrysts of euhedral olivine and plagioclase are present in some low-MgO samples (Fig. 5F). The presence of olivine and plagioclase phenocrysts in lower-Mg samples and the abundant eruptive olivine, plagioclase, and small subhedral clinopyroxene imply a crystallization sequence of olivine → plagioclase → clinopyroxene. This sequence is characteristic of the low-pressure crystallization of a low-H₂O melt.

Major elements, Cr, Ni, and V were analyzed using fused beads while Sc, Rb, Sr, Zr, Nb, and Y were analyzed using pressed pellets. All were analyzed via X-ray fluorescence (XRF) using a Philips PW2400 4kW automated XRF spectrometer system with a rhodium 60 kV end window X-ray tube. All other trace elements were analyzed via inductively coupled plasma–mass spectrometry (ICP-MS) using a Perkins Elmer quadruple ICP-MS following multi-acid dissolution. Nd and Sm isotopes were analyzed via thermal ionization mass spectrometry (TIMS) on a Thermo Triton™ mass spectrometer following multi-acid dissolution and column chromatography. Detailed descriptions of methods can be found in the GSA Data Repository.¹

Sample F624B was collected at the base of the Hula Hula diamictite (Fig. 4A) and consists of gravel-sized lithic fragments in a fine-grained greenschist-facies volcaniclastic matrix. There is no evidence for a glaciogenic origin of the volcaniclastic rocks at or below this horizon; however, this volcaniclastic unit is immediately overlain by stratified, matrix-supported diamictite with outsized, foreign clasts of carbonate and quartzite (Macdonald, 2011).

Approximately 2 kg of sample were crushed by jaw crusher and milled with a disk grinder. Accessory minerals were separated with a Wilfley table and Frantz magnetic separator and were concentrated by heavy liquid separation in methylene iodide (MEI) at Harvard University. Zircon was annealed at 900 °C for 60 h in a muffle furnace, mounted with natural zircon standards in a 2.5 cm epoxy round mount, polished, and imaged by cathodoluminescence (CL) prior to laser-ablation (LA) ICP-MS analysis, which was undertaken at the University of Arizona LaserChron Center utilizing the Nu Instruments HR ICP-MS and attached Photon Machines Analyte G2 Excimer laser following protocols of Gehrels et al. (2008). Common Pb corrections were made using ²⁰⁴Hg-corrected ²⁰⁴Pb measurements for each analysis and initial Pb compositions of Stacey and Kramers (1975). U and Th concentrations and Pb/U fractionation were calibrated against Sri Lanka (SL) zircon standard (563.5 ± 3.2 Ma; ∼518 ppm U and 68 ppm Th; Gehrels et al., 2008). Pb/U fractionation was monitored by a secondary standard as well (R33: 421 Ma; Black et al., 2004; Mattinson, 2010). Data reduction followed Gehrels et al. (2008).

Following LA-ICP-MS analysis, selected zircon grains were removed from epoxy mounts for chemical abrasion–isotope dilution–thermal ionization mass spectrometry (CA-ID-TIMS) at Boise State University. The details of CA-ID-TIMS analysis follow those previously published by Davydov et al. (2010). U-Pb dates and uncertainties were calculated from measured isotopic ratios using the algorithms of Schmitz and Schoene (2007), the U decay constants of Jaffey et al. (1971), and the EARTHTIME ET535 tracer (Condon et al., 2007). The isotope ratios, dates, and further data reduction details are compiled in Table DR3 (see ^{footnote 1}).

Weighted mean ²⁰⁶Pb/²³⁸U dates were calculated from isotopically equivalent crystals. Errors on the weighted mean dates are the internal errors based on analytical uncertainties only, including counting statistics, subtraction of tracer solution, and blank and initial common Pb subtraction. See supplementary materials for further details (see ^{footnote 1}).

Forty-seven of 80 grains from sample F624B analyzed by LA-ICP-MS yielded Cryogenian ages (from 697 to 773 Ma; Fig. 6A) with a weighted mean ²⁰⁶Pb/²³⁸U age of 719 ± 3 Ma (mean square of weighted deviates [MSWD] = 1.8). The older grains are consistent with detrital components observed in the subjacent clastic units (Macdonald et al., 2009; Strauss et al., 2013). Of the Cryogenian grains, eight grains were selected from across the Cryogenian probability density function (Fig. 6B), and analysis by CA-ID-TIMS yielded two clusters of dates at ca. 719.5 and ca. 723.6 Ma (Fig. 6C). Six of these grains were statistically equivalent with a weighted mean of 719.47 ± 0.29 Ma (MSWD = 0.57, probability of fit = 0.989), which is interpreted as the igneous crystallization age of one of the dominant contributors to the epiclastic volcanic provenance. The weighted mean date most conservatively represents a maximum age of deposition of these volcaniclastic sediments. Given the stratigraphic proximity to the Kikiktat volcanics, it is highly probable that these grains represent a silicic component of the volcanic province, and, along with the inherited 723.6 Ma component, they establish its likely age at ca. 720 Ma, coeval with the Franklin large igneous province.

Major, trace, and normative compositions are shown in Table 1, and isotopic data are found in Table 2. Whole-rock major-element chemistry for the Kikiktat volcanics reveals that these lavas range between 48% and 52% SiO₂, while MgO concentrations vary from 5% to 12%. Continuous major-element variations (Fig. 7) indicate that substantial geochemical differentiation has occurred within this volcanic suite. Typical total alkali versus silica classification (Le Bas et al., 1992) identifies the Kikiktat volcanics as basalts, while normative mineralogy reveals them to be olivine tholeiites. Increasing FeO and TiO₂ as a function of MgO is clearly apparent (Fig. 7), indicative of a tholeiitic differentiation trend, which is reinforced by Fe_4.0 / Fe_8.0 ratios (Fe_4.0 / Fe_8.0 = tholeiitic index [THI]) of 1.23 (values greater than 1 indicate tholeiitic; values less than 1 indicate calc-alkaline). Zimmer et al. (2010) established an empirical relationship between Fe_4.0/Fe_8.0 ratios and H₂O content; based on this empirical calculation, Fe_4.0/Fe_8.0 ratios equate to ∼1% H₂O, and the estimated error on this method is ±1.2%.

Trace-element concentrations on these basalts exhibit large ion lithophile element (LILE) and light rare earth element (LREE) enrichment (Fig. 8) compared to normal mid-ocean-ridge basalts (N-MORB), whereas the heavy rare earth elements (HREEs) are moderately depleted (Sun and McDonough, 1989). Within the suite, fractionation between the LREEs and HREEs is limited, with [La/Lu]_PM ratios normalized to primitive mantle composition (Sun and McDonough, 1989) ranging from 1.7 to 3.8. The HREEs are also relatively unfractionated, with normalized [Sm/Lu]_PM ratios of ∼1.5. Trace-element anomalies (Figs. 9A–9D) are present, including positive Pb anomalies, negative Nb anomalies and a progression from positive Sr anomalies to negative Sr anomalies with decreasing MgO. Isotopically, the most primitive compositions of Kikiktat volcanics have an initial ε_Nd value of +3.91, whereas low-MgO samples have values as low as -12.

Based on the International Union of Geological Sciences (IUGS) high-Mg classification scheme (Fig. 10; Le Bas, 2000), the most primitive Kikiktat samples range from picrites to picro-basalts, whereas the more-evolved rocks are basalts. This picritic affinity of the high-MgO samples is supported by the petrographic characteristics described earlier, including abundant olivine phenocrysts and spherulitic and variolitic clinopyroxene (Kerr and Arndt, 2001).

Continuous major- and trace-element variations indicate that the compositional variations record the evolution of a liquid through fractional crystallization. Plots of compatible (i.e., Ni) versus incompatible (i.e., Zr) minor elements (Fig. 11) reveal substantial decreases in compatible element abundances coinciding with increasing incompatible element abundances. This is consistent with fractional crystallization rather than a partial melting trend (Cocherie, 1986). However, the overall variations in Zr abundances (45 ppm to 170 ppm) could possibly imply unrealistic degrees of fractional crystallization. Consequently, some degree of assimilation and fractional crystallization (AFC) may be required (Fig. 11).

Within an established AFC framework coupled with the observed mineralogy, an empirical model can be established to explain the major-element compositional variations observed in the Kikiktat volcanics. Relatively constant Fe concentrations of the most primitive samples (i.e., >∼8% MgO) with decreasing MgO (Fig. 7A) suggest a dominant role for early olivine fractionation. More-evolved samples (i.e., <∼8% MgO) show increasing FeO (and TiO₂) with decreasing MgO, characteristic of a tholeiitic trend of iron enrichment. This inflection in the FeO and TiO₂ liquid line of descent (Figs. 7A and 7D at ∼8%–9% MgO) most likely coincides with the appearance of plagioclase as the second fractionating phase. Additionally, FeO and TiO₂ undergo continued enrichments across the entire spectrum of MgO compositions, indicating that Fe-Ti oxides were never liquidus phases. This logic can be applied to P variations, which indicate that apatite saturation was also never reached. Ca concentrations lie below those of typical of MORB, even at primitive compositions before the appearance of diopside and/or Ca-rich plagioclase. The most primitive melt is in equilibrium with an olivine of ∼Fo₈₈ composition (calculated using a Kd = 0.3 [Roeder and Emslie, 1970] and Fe²⁺/Fe_T = 0.85 based on V/Sc ratios [Lee et al., 2005]), has not undergone extensive AFC, and has a CaO concentration of only 9.75%. This composition plots outside the lherzolite source field (Fig. 7E) defined by Herzberg and Asimow (2008) and in the low-Ca field associated with a pyroxenite or harzburgite source (i.e., clinopyroxene-poor). This harzburgite association is also present in Al₂O₃-SiO₂ plots (Fig. 7F), with high-MgO samples plotting close to the clinopyroxene-out boundary for a lherzolite source, implying either large degrees of partial melting or melting of a low-Al source (i.e., harzburgite).

Thermodynamic modeling supports this empirical scenario of olivine fractionation followed by plagioclase then clinopyroxene crystallization. In Figures 7A–7G, we have used alphaMELTS (Ghiorso and Sack, 1995; Smith and Asimow, 2005) to model the pure fractional crystallization of the highest MgO liquid. This model crystalizes olivine first at 1300 °C, followed by plagioclase and then clinopyroxene. The total compositional range represented by the Kikiktat volcanics would require ∼60% fractional crystallization. Adjusting this model through the addition of 1% H₂O results in the significant delay of plagioclase, resulting in a modeled LLOD (liquid line of decent) that does not match the observed compositional variations (Fig. 7A), reinforcing petrographic observations for an essentially dry melt. While the variations in major-element compositions are successfully explained both empirically and thermodynamically, this pure fractional crystallization model does have problems. First, Philpotts and Carroll (1996) showed >35% fractional crystallization is difficult to accomplish due to the rapid increase in strength of the crystal mush due to plagioclase and pyroxene crystallization. Second, such a model does not provide a reasonable explanation for the isotopic or trace-element variations.

Elemental concentrations in the Kikiktat volcanics are similar to those observed in oceanic plateau basalts (Kerr, 1995), with their relatively flat unfractionated trace-element profiles (Fig. 8). Significantly, the unfractionated HREEs indicate a garnet-free source (Fig. 8) and plot within the shallow melting array of Pearce (2008), having distinctly low TiO₂/Yb ratios (Fig. 12A).

Negative Sr and Eu anomalies (Sr* and Eu*) are present and covary with Zr and MgO, indicating that their development was a function of decreasing temperature and hence crystal fractionation (Figs. 9A and 9B). Both Sr* and Eu* are statistically correlated (Pearson R = 0.75); consequently, plagioclase fractionation was the likely driver of these anomalies. Negative Nb (Nb*) and positive Pb (Pb*) anomalies also occur (Figs. 9C and 9D), but unlike Eu* and Sr*, they do not covary with any proxies for fractionation (i.e., Zr and/or MgO). The sizes of the Nb and Pb anomalies are essentially random, with the only minor exception being that Pb* is consistently positive, while Nb* is consistently negative.

From an isotopic perspective, the most primitive composition of the Kikiktat volcanics has an initial ε_Nd value of +3.91, compared to an estimated depleted mantle (DMM) value of ∼+8.35 at ca. 720 Ma (Goldstein et al., 1984). Considerable variation within the suite exists, with samples as low as −12 (Table 2). Such variation in initial ε_Nd compositions again implies either crustal assimilation and/or source metasomatism during the petrogenesis of these basalts. It is important to note that the relationship between MgO, Zr, or other measures of fractionation and ε_Nd is minimal (Fig. 13A; Table DR1 [see ^{footnote 1}]), except for a weak relationship with LREE/HREE ratio (Fig. 13B). It has been argued that the heat of fusion should lead to assimilation and fractional crystallization (AFC) being strongly coupled (DePaolo, 1981); however, with respect to the Kikiktat volcanics, it is clear that AFC is not strongly coupled (Fig. 13A), with cross plots between ε_Nd and Zr being essentially random. Cribb and Barton (1996) showed that decoupled AFC processes produce differing results than coupled AFC, specifically, higher compatible element abundances; however, the observed trace-element trends are still coherent. Considering the highly variable nature of ε_Nd as a function of Zr, the contaminant may be multicomponent in nature.

Nd model ages are displaced from their ca. 720 Ma emplacement age and range between 2.32 and 3.88 Ga (average = 3.16 Ga; 1σ = 0.48). Numerous possible scenarios can be envisaged to explain these model ages, including but not limited to multiple episodes of reworking, an isotopically evolved mantle source, or material representing a mixing age between ca. 720 Ma DMM melt sources and Paleoarchean crustal contaminants (Figs. 13C and 13D; Arndt and Goldstein, 1987).

As discussed already, the variations in initial ε_Nd values, Nd model ages, Pb*, and Nb* are unlikely to be the result of pure fractional crystallization but instead may require some degree of crustal assimilation, complicated by the likelihood that the contaminant is possibly multicomponent in nature. Consequently, it is important to ascertain what elements are least affected by contamination.

As a first-order test on the effect of assimilation on trace- and major-element compositions, we performed a principal component analysis (PCA) of the data. The interpretative value of PCA over correlation coefficients, least squared regression analysis, and cross-plots is that relationships between the principal components and the variables can be quantified (i.e., the contribution of the variables to the principal components); consequently, the relative effect of assimilation on individual elemental compositions can be quantified. Utilizing the major element variables of SiO₂, TiO₂, Al₂O₃, FeO(T), MgO, CaO, Σalkalis, ΣHFSEs (high field strength elements), ΣLREEs, ΣHREEs, LREE/HREE ratios, and the isotopic variable of ε_Nd, and assuming that variations in ε_Nd are principally governed by assimilation, we statistically determined the relationship between assimilation and elemental compositions.

This analysis suggests that ∼72% of the major- and trace-element variability can be explained by just two principal components (Fig. 14). Principal component 1 (PC1), which accounts for ∼49% of the variability in the data set, is dominated by the bulk major-element compositions in which MgO and CaO show a strong positive relationship with each other but an equally strong inverse relationship with all other oxide and trace-element compositions. Considering that fractional crystallization in tholeiitic suites results in decreasing MgO and CaO compositions, while other major oxides and incompatible elements increase, we conclude that PC1 is recording mineral fractionation. Significantly, and as would be predicted if PC1 is indeed recording mineral fractionation and if assimilation and fractional crystallization are largely decoupled processes, the contribution of ε_Nd to PC1 is minimal (∼0.6% contribution). Principal component 2 (PC2), which accounts for ∼22% of the total variability in elemental and isotopic compositions, is dominated by variations in ε_Nd (25% contribution), changes in LREE/HREE ratios (22% contribution), and changes in LILEs (14% contribution), with lesser contributions from SiO₂, Al₂O₃, and total alkalis (24% total combined contribution). A first-order inference from this analysis is that while most major elements have been affected by assimilation, compositional variations are still dominantly a function of mineral fractionation (Fig. 14), which is consistent with the coherent liquid line of descent defined by the major-element data (Fig. 7). It appears as though LILEs, LREE/HREE ratios, SiO₂, Al₂O₃, and total alkalis are most affected by assimilation, which is not surprising if we consider that the likely contaminant is country rock with a similar composition to the upper continental crust (Taylor and McLennan, 1995) or metasedimentary rocks (Nance and Taylor, 1976; Taylor and McLennan, 1985). A list of the contributions of each variable to principal components 1 and 2 can be found in Table 3.

A significant result of this analysis is that LREE/HREE ratios and LILE concentrations show an overwhelming negative correlation with ε_Nd. This implies that increasingly evolved ε_Nd values are associated with increasing LREE/HREE ratios (Fig. 13B) and increasing overall LILE abundances. Therefore, the variations in ε_Nd, LREE/HREE ratios, and LILE abundances are most likely cogenetic, and the contaminant is LILE- and LREE-enriched. Hence, assimilation, while a significant factor in the petrogenesis of this suite of lavas, has not entirely obscured the primary major- (i.e., MgO, CaO, TiO₂) and trace-element (i.e., HFSE and HREE) characteristics of the primitive high-MgO samples. This is an important consideration in applying model calculations of the parental melt to inferring primary melt compositions, mantle potential temperatures, and source composition.

Based on our assessment of the effects of assimilation on major-element compositions, high field strength element ratios and abundances, and the recognition that the high-MgO samples have only fractioned olivine, we can use phase equilibrium constraints for the picritic parental melt to investigate a potential mantle source for the Kikiktat volcanics.

Using alphaMELTS (Ghiorso and Sack, 1995; Smith and Asimow, 2005), we modeled the phase diagram for the parental melt (Fig. 15A). Based on the observed crystallization sequence of olivine, plagioclase, and clinopyroxene, we surmise that magma storage, subsequent fractional crystallization, and assimilation processes occurred at crustal levels of less than ∼1.6 kbar or ∼5 km depth. Considering this shallow depth, the primitive nature of this parental melt, and the consensus that the associated sedimentary succession records Kikiktat-age rifting (Clough and Goldhammer, 2000; Macdonald et al., 2009; Strauss et al., 2013), taken together, these observations provide evidence in support of significant crustal attenuation by this time. While many voluminous basaltic occurrences have been attributed to a deep mantle melting in the absence of a garnet signature, we would argue that the most parsimonious explanation for the absence of a garnet signature in the rare earth element profiles for these lavas (Figs. 8 and 12) is a result of mantle melting at depths outside of garnet stability (i.e., less than ∼27 kbar [∼81 km]). In isolation, this pressure limit may remove any scenario that involves a “pure” plume source for these melts, irrespective of calculated mantle potential temperatures.

Arguably the most significant aspect of the phase equilibrium constraints for the parental melt is that clinopyroxene is never in equilibrium with the parental liquid at any pressure. While this conclusion is based upon a parental melt in equilibrium with olivine of Fo₈₈ composition, these conclusions remain the same if we use estimates for the composition of the primary melt. For the moment, ignoring dynamic melting processes (e.g., Eggins, 1992), these phase constraints result in the parental melt being in equilibrium with harzburgite at pressures of ∼3.2 kbar or depths of ∼10 km. In comparison, typical basalts (i.e., primitive MORB) have olivine, clinopyroxene, and orthopyroxene as liquidus phases (Green and Ringwood, 1967; Pertermann and Hirschmann, 2003; Yasuda et al., 1994). Unlike the MORB comparison, this harzburgite relationship is shared with primitive melts associated with the contemporary Natkusiak lavas (Fig. 15B) of Victoria Island (Dostal et al., 1986). Considering their practically identical ages and hypothesized paleogeographic proximity with restoration somewhere along the northern margin of Laurentia (see discussion in Macdonald et al., 2009; Strauss et al., 2013), their petrogenetic similarities naturally support a link between these volcanic occurrences as both being part of the Franklin large igneous province. However, it must be noted that the Franklin large igneous province covers an area of over 2.25 × 10⁶ km² (Ernst et al., 2008), with lavas, sills, and dikes extending over much of northern Canada from Yukon to Greenland (Denyszyn et al., 2009; Dostal et al., 1986), such that more precise links to a particular province of the Franklin large igneous province require additional constraints.

The volume and spatial extent of the Natkusiak lavas is far greater (Dostal et al., 1986; Dupuy et al., 1995; Shellnutt et al., 2004) than the Kikiktat volcanics. Their present-day thickness is ∼1100 m (Dostal et al., 1986), representing the thickest extrusive equivalent of the Franklin large igneous province preserved along the Arctic margin. Dostal et al. (1986) argued for a lherzolite source for the Natkusiak lavas, whereas Shellnutt et al. (2004) and Dupuy et al. (1995) also suggested an upper-mantle source (i.e., lherzolite). Consequently, the phase relationships for these lavas (Fig. 15B) require either extensive melting of a highly fertile lherzolite source to leave behind a harzburgitic residue free of clinopyroxene, or they represent lessor degrees of partial melting of a more-depleted source. Considering the voluminous nature of the Franklin lavas, sills, and dikes, a relatively fertile lherzolite source may be favored.

Accepting the conclusions that the Natkusiak lavas were derived from a lherzolite source (Dostal et al., 1986; Dupuy et al., 1992; Shellnutt et al., 2004), geochemical comparisons between the Natkusiak lavas and Kikiktat volcanics can be used to constrain the nature of the source and percentage of melt of the latter. To this end, we apply the Mg8.0 method (Klein and Langmuir, 1987; Turner and Hawkesworth, 1995), which involves the recalculation of major-element data to the common MgO content of 8%. The goal of this method is to be able compare disparate lava compositions at a similar point in their petrogenesis; since MgO has a strong linear relationship to temperature (e.g., Niu et al., 2002), MgO is the basis of the regression. While the choice of 8% MgO is a somewhat arbitrary comparison point, it does broadly coincide with liquids hitting the olivine-plagioclase-clinopyroxene cotectic; consequently, two linear regressions (either side of 8% MgO) can be undertaken to adjust compositions back to 8% MgO. Importantly, compositions are recalculated using equations of curves fitted to the data (typically least squared regressions), and so this type of analysis neither assumes nor requires the phase proportions associated with the crystallizing assemblage.

Figure 16 demonstrates that at 8% MgO, the major-element compositions of the Natkusiak and Kikiktat volcanics are distinct. Al_8.0 and Si_8.0 compositions for the Kikiktat volcanics are displaced to higher Al and lower Si values, consistent with a smaller degree of melting with respect to the Natkusiak lavas. Furthermore, Fe_8.0 and Ti_8.0 support a smaller degree of melting, but they also imply that the source for the Kikiktat volcanics was more depleted than that of the Natkusiak lavas. Na_8.0 contradicts this apparent source depletion; however, considerable complexity exists in interpreting Na_8.0 when the pressure of melting and/or percentage melting is not comparable. While depleted sources should produce melts with lower Na content, this effect can be masked by the increased Na concentrations associated with small percentage melts. Another factor that can complicate Na content is the significant changes in the incompatibility of Na in clinopyroxene as a function of pressure (Langmuir et al., 1992). All of this is further complicated by Na mobility under greenschist-facies metamorphism.

The compositional distinctness of the Kikiktat volcanics, with respect to low Ti and Fe and high Si (characteristics of a depleted source) and relatively high Al and Na (with respect to the Natkusiak lavas), implies that these volcanics possibly represent small percentage melting of a harzburgitic source. This interpretation is supported by high ¹⁴⁷Sm/¹⁴⁴Nd ratios, where progressive melting of mantle peridotite from lherzolite to harzburgite will be accompanied by increasing ¹⁴⁷Sm/¹⁴⁴Nd ratios of the residual source, and basalts derived from such a depleted source will consequently inherit elevated ¹⁴⁷Sm/¹⁴⁴Nd ratios. The most primitive ¹⁴⁷Sm/¹⁴⁴Nd ratios of the Kikiktat volcanics (Fig. 13C) are higher than most MORB and oceanic-island basalt examples and significantly higher than the majority of continental flood basalts (Goldstein et al., 1984)—this is despite the fact that assimilation drives the melt compositions to lower ¹⁴⁷Sm/¹⁴⁴Nd ratios (Goldstein et al., 1984).

Considering that the parental melts for neither the Natkusiak basalts nor the Kikiktat volcanics contained clinopyroxene as a liquidus phase at any pressure, and that the former required a lherzolite source (Dostal et al., 1986), the Natkusiak basalts must represent melting past the clinopyroxene-out boundary. In contrast, the Kikiktat volcanics must represent a smaller percentage melt of a harzburgitic source. Consequently, the significant compositional contrasts between the two can be reconciled with their broadly comparable phase equilibria, which predict liquids that do not coexist with clinopyroxene.

Estimating a primary melt composition for these rocks is difficult because olivine phenocrysts exhibit significant iddingsite alteration. Consequently, correcting the melt back to the primary composition of coexisting olivine is problematic. Nevertheless, calculating a primary melt and associated mantle potential temperature (T_p) is informative because mantle T_p calculations can tie these volcanics to larger regional magmatic events and support or invalidate a model for a harzburgite source. As discussed previously, the petrology and geochemistry of the Kikiktat volcanics imply relatively dry melting; along with constraints from harzburgite melting experiments (Falloon and Danyushevsky, 2000; Falloon et al., 1988), this indicates a higher mantle T_p for the generation of these melts than conventionally recognized for the upper mantle.

Two methods have been applied to achieve this primary melt calculation. First, we applied the PRIMELT2 model of Herzberg and Asimow (2008), which couples an inverse model (i.e., olivine addition) with a forward melting model, such that the intersection of the two indicates the primary mantle melt. Considering the difficulty highlighted above in regard to primary olivine compositions, the PRIMELT2 method does not assume nor require a forsterite composition of coexisting olivine; furthermore, the method highlights compositions that may be a result of pyroxenite melting and/or melting in the presence of CO₂ (Herzberg and Asimow, 2008). However, the PRIMLET2 method does neglect the effects of assimilation, possibly introducing errors into the T_p calculations. The second method used is an extension of the Mg_8.0 method, where only those samples with greater 9% MgO are used to calculate a least-squared regression line. The regression is then extrapolated to higher MgO values until mantle olivine compositions are reached, essentially producing the reverse of assimilation and fractional crystallization.

In Figure 17, the PRIMELT2 method estimates a picritic primary melt with MgO = 15.6%, which requires ∼11% olivine addition and results in equilibrium olivine compositions of Fo_91.2 indicative of a highly depleted peridotite source. Calculated mantle potential temperatures of 1475 °C are obtained via the PRIMELT2 method. Also shown in Figure 17 is the reverse AFC model, which, if extended back to the composition of Fo_91.2 calculated using PRIMELT2, yields a calculated MgO of 15.7% for the primary melt and a mantle potential temperature of 1476 °C. This result highlights that fact that the small amount of assimilation in this model has very little effect on the primary melt calculation and associated mantle potential temperatures. This conclusion is not surprising when considering that mantle potential temperatures are primarily a function of MgO (Herzberg and Asimow, 2008; Herzberg et al., 2007; Iwamori et al., 1995) and that the reverse AFC vector and olivine addition vector (Fig. 17) diverge by only a small amount in MgO space.

One of the distinguishing features of the Kikiktat volcanics is their low TiO₂ content. This unique composition is not a result of fractionation, as Fe-Ti oxides are never a liquidus phase, and tholeiitic fractionation trends drive TiO₂ compositions to higher not lower concentrations. While assimilation of upper continental crust can result in lower TiO₂ compositions, it can be seen in Figures 7 and 14 that TiO₂ compositions are not overwhelmingly affected by assimilation, especially the more primitive samples (samples with >9% MgO). In fact, the TiO₂ contents of the primitive melts are lower than N-MORB compositions, and neither the primitive nor the evolved melts record the significant negative Ti anomalies (Ti*) associated with assimilation of upper continental crust (Fig. 8; Sun and McDonough, 1989; Taylor and McLennan, 1995). Consequently, the TiO₂ of the pre-AFC melt must be a function of source composition and contain unusually low Ti abundances.

While these low-Ti Kikiktat volcanics are presumably associated with the breakup of Rodinia, continental flood basalts associated with the breakup of Pangea also have spatially distinct low-Ti magmatic provinces. For example, the Gondwanan Ferrar low-Ti province stretches some 3000–4000 km across Australia and Antarctica (Hergt et al., 1991). Furthermore, the Parana continental flood basalts of Brazil are an example of a smaller, but still significant Gondwanan low-Ti province (Gibson et al., 1995). Gibson et al. (1995) argued that the combination of low Ti, Fe, and P and an evolved ε_Nd signature requires a mantle source for the Ferrar that was depleted in major and trace elements, with the subcontinental lithospheric mantle being the most likely candidate. In contrast, the combination of coexisting low-Ti and high-Ti basalts for the Parana continental flood basalt indicates plume–subcontinental lithospheric mantle interaction (Gibson et al., 1995). Following similar arguments, the Kikiktat volcanics can be plausibly linked to melting of subcontinental lithospheric mantle.

An interesting question that follows is: Was this depleted subcontinental lithospheric mantle a preexisting feature of the Laurentian lithosphere, or is it a consequence of the mantle plume suggested for the Franklin large igneous province? The thermochemical plume model (Campbell, 2007; Campbell and Griffiths, 1990; Griffiths and Campbell, 1990; Watson and McKenzie, 1991; Wyllie, 1988) suggests that hot, adiabatically decompressing mantle (plume head) impinges on and flattens at the base of the lithosphere, and so most melting should occur in the plume head, resulting in flood basalt volcanism. Such extensive melting associated with this flood magmatism would consequently leave a melt-depleted (and lower--density) mantle that may flow laterally away from the plume center. Consequently, volcanism associated with the Kikiktat volcanics may represent second-stage melting of the Franklin plume head at the “far-flung” fringes of the magmatic province. However, considering the debate that surrounds the natural variability in mantle potential temperatures (Dalton et al., 2014; Herzberg, 2011; Herzberg et al., 2007), it is safer to assert that the precise tectonic-thermal regime for the Kikiktat volcanics remains ambiguous.

During the breakup of Gondwana, the emplacement of continental flood basalts through crust with different subcontinental lithospheric mantle histories resulted in magmas with different compositions (Gibson et al., 1995). Specifically, low-Ti basalts are associated with subcontinental lithospheric mantle regions further afield from the plume head, whereas high-Ti basalts are associated with areas that experienced melting directly associated with the plume head. Applied to the Franklin large igneous province, this would suggest that the high-Ti basalts of Victoria Island are possibly associated with eruption close to the plume head, whereas the Kikiktat volcanics formed further afield in a low-Ti province. Based on regional stratigraphic evidence for doming prior to the emplacement of the Natkusiak basalts (Rainbird, 1993) and the orientation of radiating dike swarms (Ernst and Bleeker, 2010), the Franklin plume head has been tied to a location near the northwestern margin of Laurentia in the vicinity of northern Victoria Island (Fig. 18A).

Ignoring pre-Mesozoic displacement of terranes along the Arctic margin, the simple counterclockwise rotation model for the opening of the Arctic Ocean (review by Gottlieb et al., 2014; Lawver and Scotese, 1990) predicts that (Macdonald et al., 2012) the North Slope of Arctic Alaska would restore off-board of Banks and Victoria Islands. If this simple restoration is correct, the expectation is that the Kikiktat volcanics would share similar petrogenetic and geochemical relationships with the high-Ti Natkusiak basalts. However, this is inconsistent with the data presented herein and can be reconciled either with strike-slip displacement during the Mesozoic opening of the Arctic Ocean or by pre-Mesozoic movement of the North Slope subterrane along the Arctic margin. A Neoproterozoic restoration (Fig. 18B) of the North Slope subterrane between Siberia and Laurentia (Ernst et al., 2012; Macdonald et al., 2009) along the northern margin of Greenland or Ellesmere Island and further away from the Franklin plume head has also been suggested by stratigraphic, geochronological, and paleontological data (Strauss et al., 2013).

The 719.5 ± 0.3 Ma maximum depositional age obtained from volcaniclastic strata directly below the Hula Hula diamictite is consistent with both a Sturtian age for the Hula Hula diamictite and previous age constraints on the onset of the Sturtian glaciation in Yukon, Canada, between 717.4 ± 0.1 and 716.5 ± 0.2 Ma (Macdonald et al., 2010) and in China around 715.9 ± 2.8 Ma (Lan et al., 2014). Assuming the abundant epiclastic minerals originated from proximal felsic lavas or tuffs of the underlying Kikiktat volcanics, these 723.5–719.5 Ma zircon grains provide the best estimate for the age of the volcanic province. The confirmation of the intimate relationships among the breakup of Rodinia, the emplacement of the Franklin large igneous province, and the onset of the Sturtian glaciation implies a mechanistic relationship. This is consistent with recent work that suggests that the weathering of widespread magmatic provinces associated with the polyphase breakup of Rodinia could have increased global weatherability, as seen in Sr, Os, and Nd isotope systematics (Halverson et al., 2010; Rooney et al., 2014), leading to effective CO₂ consumption and the initiation of global glaciation. As more occurrences of the Franklin large igneous province are identified, both in peri-Laurentian terranes such as Arctic Alaska and on other margins such as Siberia (Ernst et al., 2012), further estimates will be needed to detail the full extent of the Franklin large igneous province and its effects on climate and ocean geochemical cycles.

The Kikiktat volcanics represent a suite of well-preserved continental flood basalts, including picrites and picro-basalts, that have undergone shallow tholeiitic fractionation from a picritic parental magma. While assimilation has played a role in their petrogenesis, the effect is not overwhelming. The combination of low Ti, Ca, and Fe with moderate Al and Na compositions can best be explained by derivation from a harzburgitic mantle involving the subcontinental lithospheric mantle. Considering the elevated model mantle T_p (1476 °C), low-Ti characteristics, and their indistinguishable ages, we suggest that the Kikiktat volcanics represent a low-Ti province of the much larger Franklin large igneous province. The new 719.5 ± 0.3 Ma CA-ID-TIMS age on volcaniclastic rocks at the top of the Kikiktat volcanics presented herein provides a new constraint on the timing and extent of the Franklin large igneous province and an additional maximum age constraint on the onset of the Sturtian glaciation, further linking these events in time and space.

This research was supported by Natural Sciences and Engineering Research Council of Canada (NSERC) Discovery grants to Halverson and a NSERC Vanier Fellowship to Cox. Strauss was supported by a National Science Foundation (NSF) Graduate Research Fellowship; W.C. McClelland acknowledges support from NSF Tectonics (EAR-1049368 and EAR-1032156 to the University of Arizona LaserChron Center); and Macdonald acknowledges support from NSF Tectonics (EAR-1049463), NSF Sedimentary Geology and Paleontology (EAR-1148058), and the Massachusetts Institute of Technology NASA Astrobiology Institute node. We would like to thank Charlie Langmuir and Zhongxing Chen for use of the inductively coupled plasma–mass spectrometry laboratory at Harvard University; Mark Pecha and Blake Hodgin for analytical support; Tim Gibson, Esther Kennedy, and Alex Gould for field assistance; and Kirk Sweetsir of Yukon Air Service, Christian Cabanilla of ERA Helicopters, and Nancy Brandt and Dana Troufler-Moudra of Polar Field Services for safe transportation within the Brooks Range. Anne Marie de la Rosa of the U.S. Fish and Wildlife Service enabled our access to the Arctic National Wildlife Refuge.