Accurate reconstruction of the Late Cretaceous paleogeography and tectonic evolution of the western North American Cordilleran margin is required to resolve the long-standing debate over proposed large-scale, orogen-parallel terrane translation. The Nanaimo Basin (British Columbia, Canada) contains a high-fidelity record of orogenic exhumation and basin subsidence in the southwestern Canadian Cordillera that constrains the tectonic evolution of the region. Integration of detrital zircon U-Pb geochronology, conglomerate clast U-Pb geochronology, detrital muscovite 40Ar/39Ar thermochronology, and Lu-Hf isotopic analysis of detrital zircon defines a multidisciplinary provenance signature that provides a definitive linkage with sediment source regions north of the Sierra Nevada arc system (western United States).

Analysis of spatial and temporal provenance variations within Nanaimo Group strata documents a bimodal sediment supply with a local source derived from the adjacent magmatic arc in the southern Coast Mountains batholith and an extra-regional source from the Mesoproterozoic Belt Supergroup and the Late Cretaceous Atlanta lobe of the Idaho batholith. Particularly robust linkages include: (1) juvenile (εHf >+10) Late Cretaceous zircon derived from the southern Coast Mountains batholith; (2) a bimodal Proterozoic detrital zircon signature consistent with derivation from Belt Supergroup (1700–1720 Ma) and ca. 1380 Ma plutonic rocks intruding the Lemhi subbasin of central Idaho (northwestern United States); (3) quartzite clasts that are statistical matches for Mesoproterozoic and Cambrian strata in Montana and Idaho (northwestern United States) and southern British Columbia; and (4) syndepositional evolved (εHf >−10) Late Cretaceous zircon and muscovite derived from the Atlanta lobe of the Idaho batholith. These provenance constraints support a tectonic restoration of the Nanaimo Basin, the southern Coast Mountains batholith, and Wrangellia to a position outboard of the Idaho batholith in Late Cretaceous time, consistent with proposed minimal- fault- offset models (<~1000 km).

Sedimentologic and stratigraphic analyses of Mesozoic forearc basins along the western Laurentian continental margin provide outstanding first-order constraints on the tectonic evolution of the margin (e.g., Dickinson, 1976, 1995; Dickinson and Seely, 1979; Ingersoll, 1979, 1983; Busby-Spera and Boles, 1986; Mustard, 1994). Regional analysis of stratigraphic successions within these forearc sequences relies on litho- and biostratigraphic correlations, which are useful but limited in scope. Detrital zircon geochronology provides a complementary method of basin analysis that permits direct linkages between basinal strata and their source regions (e.g., Gehrels et al., 1995; Gehrels, 2000; Kimbrough et al., 2001; DeGraaff-Surpless et al., 2003; Fedo et al., 2003; Surpless et al., 2006; Trop and Ridgway, 2007). Increasingly sophisticated provenance tools, such as large-n detrital zircon geochronology (Pullen et al., 2014; Matthews and Guest, 2016), detrital thermochronology (Reiners et al., 2005; Reiners and Brandon, 2006; Carrapa et al., 2004), and hafnium isotopic analysis (Bahlburg et al., 2009; Gehrels and Pecha, 2014), permit high-resolution analysis of sediment provenance and basin evolution. There are, however, numerous unresolved issues regarding the spatial and temporal evolution of Mesozoic forearc basins along the North American western margin, including controversies about the timing (Surpless et al., 2006; Haggart et al., 2006, and references therein; Ward et al., 2012; Kent et al., 2020), provenance (Matthews et al., 2017; Surpless and Gulliver, 2018; Orme and Surpless, 2019; Isava et al., 2021), paleogeographic setting (Garver and Davidson, 2015; Mahoney et al., 2016; Sauer et al., 2017a; Stevens Goddard et al., 2018; Trop et al., 2020), and tectonic evolution of different basins along the margin (Wright and Wyld, 2007; Sauer et al., 2017a, 2019; Coutts et al., 2020).

This investigation addresses the long-standing controversy between geological and paleomagnetic evidence for large-scale, orogen-parallel terrane displacements in the Canadian Cordillera (Beck and Noson, 1972). The paleogeographic setting of the Late Cretaceous Nanaimo Basin of the southwestern Canadian Cordillera has been the subject of debate for decades (Cowan et al., 1997; Mahoney et al., 1999; Matthews et al., 2017; Sauer et al., 2017a). Paleomagnetic data and provenance studies have been used to argue both for (Ward et al., 1997; Enkin, 2006; Krijgsman and Tauxe, 2006; Matthews et al., 2017; Sauer et al., 2017a) and against (Mahoney et al., 1999, 2016; Butler et al., 2001; Stamatakos et al., 2001; Kim and Kodama, 2004) large-scale latitudinal translation of the basin along the western continental margin.

Accurate reconstruction of Late Cretaceous paleogeography and tectonic evolution of the western Cordilleran margin requires a full and integrated analysis of spatial and temporal variations in basin subsidence patterns and sedimentary provenance along the margin. In this contribution and in a companion paper (Isava et al., 2021), we present an integrated model for the evolution of the Nanaimo Basin that utilizes a large number of independent data sets—sedimentology, paleontology, geochronology, thermochronology, and isotopic geochemistry—and provides critical constraints on Nanaimo Basin evolution, plate kinematics, and the dynamic linkages between plate subduction, orogenic exhumation, and basin evolution. Several independent provenance indicators, including detrital zircon U-Pb geochronology and Hf isotopic geochemistry, conglomerate clast U-Pb geochronology, and muscovite Ar-Ar thermochronology, are applied herein to establish a definitive link between the Nanaimo Basin and source regions north of the Sierra Nevada magmatic arc system (western United States). We demonstrate that these results are compatible with known Cretaceous–Paleocene fault offsets along the continental margin but at odds with proposed paleomagnetic reconstructions that indicate large-scale terrane translation. We further test the large-scale translation hypothesis by comparing proposed source regions both north and south of the Sierra Nevada arc segment and demonstrate that southern source regions are incompatible with data from the Nanaimo Group.

The western Late Cretaceous Laurentian margin was a complex collage of accreted terranes, volcanic arcs, accretionary prisms, fold-and-thrust belts, and marginal basins, all variously disrupted by widespread magmatism and structural deformation (Fig. 1; Coney, 1972; Monger et al., 1982; Dickinson, 2004). The rate and style of orogenic exhumation and basin subsidence along the margin varied widely as a function of plate kinematics and subduction geometry (Rea and Dixon, 1983; Engebretson et al., 1985; Kimbrough et al., 2001; Trop, 2008; Jacobson et al., 2011; Ingersoll, 2012; Yokelson et al., 2015; Stevens Goddard et al., 2018). Margin-parallel variations in the rate and vector of convergence during the Late Cretaceous led to latitudinal variations in tectonic underplating, lithospheric delamination, and crustal thickening that controlled magmatic episodicity and patterns of structural deformation along the Cordilleran margin (Ducea, 2001; Saleeby, 2003; Paterson et al., 2004; Jacobson et al., 2007; Gehrels et al., 2009; Mahoney et al., 2009; Miller et al., 2016).

The Canadian segment of this tectonically complex margin is a mosaic of allochthonous tectonostratigraphic terranes, juxtaposed along regional fault systems and intruded by Jurassic to Cenozoic plutons (Figs. 1, 2; Coney et al., 1980; Wheeler et al., 1991; Nelson et al., 2013; Monger and Gibson, 2019). These tectonostratigraphic terranes have traditionally been grouped into the Insular and Intermontane superterranes, interpreted to be characterized by separate and distinct geologic histories (Monger et al., 1982). Subsequent geologic studies incorporating isotopic and paleontologic data indicate a more complex geologic evolution and require revision of this model (Monger and Gibson, 2019). It is still useful to group these terranes geographically into those composing the inner, or Intermontane, portion of the Cordillera and those forming the outer, or Insular, portion of the Cordillera located along the Pacific margin (Fig. 1). To the east, the Intermontane assemblage of terranes consists of the Stikine, Cache Creek, and Quesnellia terranes, which are thought to have been accreted to the western margin of North America by Middle Jurassic time (Ghosh and Lambert, 1995; Nelson et al., 2013); to the west, the Insular assemblage comprises the Wrangellia and Alexander terranes, which have been interpreted to have been accreted to the western margin as late as Late Cretaceous time (Sigloch and Mihalynuk, 2017; Trop et al., 2020) but potentially as early as Middle Jurassic time (Fig. 1; Rubin et al., 1990; Thompson et al., 1991; van der Heyden, 1992; Nelson et al., 2013; Monger and Gibson, 2019). Several small, disparate terranes, including the Cadwallader, Bridge River, and Methow terranes, represent Paleozoic to Mesozoic basinal assemblages floored by oceanic crust that are structurally imbricated between the Intermontane and Insular assemblages (Fig. 1; Haugerud et al., 1996; DeGraaff-Surpless et al., 2003; Monger and Brown, 2016).

A major zone of mid-Cretaceous to Cenozoic deformation and plutonism, extending from Washington State (USA) to the Yukon (Canada), lies between the Intermontane and Insular assemblages and is coincident with the morphogeologic Coast Belt of western Canada (Figs. 1, 2; Wheeler et al., 1991). This zone of mid-Cretaceous to Eocene (ca. 150–45 Ma) synorogenic plutons and metamorphic rocks, flanked on both sides by deformed sedimentary and volcanic strata, has been referred to as the Coast-Cascade orogen (Monger and Brown, 2016). The orogen is structurally complex, characterized by east- and west-vergent folds, thrust faults, and reverse faults produced by orogen-normal compression that are in part coeval with orogen-parallel strike-slip faults, shear zones, and mineral fabrics in plutonic and metamorphic rocks (Hollister and Andronicos, 1997; Rusmore and Woodsworth, 1991; Angen et al., 2014). These structures are interpreted to record Middle Jurassic to mid-Cretaceous sinistral transpression followed by Late Cretaceous to Eocene dextral transpression (Monger et al., 1994; Nelson et al., 2013; Monger and Brown, 2016).

The timing and magnitude of movement of terranes and terrane assemblages along the western Laurentian continental margin during the Late Cretaceous to early Cenozoic is a fundamental controversy in Cordilleran tectonics (e.g., Cowan et al., 1997). A major discrepancy exists between paleogeographic models derived primarily from paleomagnetic data and those based on geologic data. Paleomagnetic data indicate large-magnitude translation of both the Intermontane and Insular assemblages along the continental margin between ca. 105 and 55 Ma (Enkin, 2006; Monger, 2014). Paleomagnetic studies of sedimentary and bedded volcanic strata from the Intermontane assemblage suggest northward translation of <1100 km (Irving et al., 1995; Johnston et al., 1996; Haskin et al., 2003; Enkin et al., 2006), whereas some bedded strata from the Insular assemblage indicate >3000 km of northward displacement (Wynne et al., 1995; Ward et al., 1997; Enkin, 2006). The “Baja B.C.” (where “Baja” refers to Baja California, Mexico, and “B.C.” refers to British Columbia, Canada) tectonic model was developed to explain the apparent differential motion between the Intermontane and Insular assemblages required by the paleomagnetic data (Fig. 3, top panel; Irving, 1985; Umhoefer, 1987; Irving et al., 1996; Cowan et al., 1997). The Baja B.C. model posits that the Intermontane and Insular assemblages were separate entities during mid- to Late Cretaceous time, with the Intermontane assemblage situated ∼1100 km to the south of its present location (north of the Klamath Mountains, California and Oregon, USA) in late Early Cretaceous time, and the Insular assemblage located ∼3000 km to the south (at the latitude of Baja California, Mexico) in Late Cretaceous time (Figs. 1 and 3, top panel). The model requires that both superterranes were subsequently translated rapidly northward during the Late Cretaceous, with final docking at their present latitude in the Late Cretaceous or early Cenozoic.

The Baja B.C. model is difficult to reconcile with geologic and paleontologic constraints. The lack of a readily observable crustal-scale structure capable of accommodating 1000–3000 km of latitudinal dextral displacement, either between the Intermontane and Insular assemblages or east of the Intermontane assemblage, is problematic (e.g., Monger, 2014; McMechan et al., 2021). Proposed terrane linkages between the Intermontane and Insular assemblages as far back as the Middle Jurassic contradict models involving paleogeographically distinct superterrane evolution (van der Heyden, 1992; Rusmore et al., 1988; Monger, 2014; Monger and Gibson, 2019). In addition, the biogeographic distribution of numerous and distinct faunal groups, including bivalves, ammonites, radiolarians, foraminifers, and calcareous nannofossils, suggests that the Insular assemblage had become established in a relatively high-latitude position by Middle Jurassic time (Sliter, 1972; Jeletzky, 1984; Taylor et al., 1984; Haggart et al., 1993; Haggart and Carter, 1994; Carter and Haggart, 2006; Schröder-Adams and Haggart, 2006), although macroplant remains have been interpreted to suggest both low- and high-latitude settings, depending on the plant group assessed (Pearson and Hebda, 2006; Miller et al., 2006). Detrital zircon provenance data from the Nanaimo Group have also been used to argue both for (Housen and Beck, 1999; Matthews et al., 2017; Sauer et al., 2017a, 2019) and against (Mahoney et al., 1999, 2016; Dumitru et al., 2016) large-scale latitudinal translation.

Analysis of geologic and paleomagnetic data from Upper Cretaceous strata along the Intermontane-Insular boundary in south-central British Columbia has significantly compounded the debate (Fig. 3, bottom panel). Replicated paleomagnetic studies along the boundary verify earlier estimates of large-scale latitudinal displacement for both the Insular and Intermontane assemblages, making them difficult to discount (Irving and Thorkelson, 1990; Wynne et al., 1995; Irving et al., 1996; Haskin et al., 2003; Enkin et al., 2003; Enkin, 2006). However, integrated stratigraphic, geochemical, geochronologic, and paleomagnetic investigations have documented that Upper Cretaceous strata overlap both the Insular and Intermontane assemblages by 95 Ma, thus requiring that the superterranes constituted an enormous crustal block that underwent rapid, large-magnitude dextral and sinistral translations along the continental margin during Late Cretaceous–Cenozoic time (Fig. 3, bottom panel; Enkin, 2006). Kent and Irving (2010) developed a new apparent polar wander path for the Triassic–Paleogene time frame using inclination-corrected paleomagnetic sites that supports large-scale sinistral then dextral translation of the combined Intermontane-Insular crustal block during this time frame.

Adherence to the paleomagnetic model thus requires large-scale episodic translation of the combined Intermontane-Insular crustal block, with deposition of upper Lower Cretaceous volcanic rocks (ca. 105–100 Ma; Spences Bridge Group; Haskin et al., 2003) at the latitude of Oregon followed by sinistral translation of this enormous crustal block ∼2000 km to the south to the latitude of Baja California (Fig. 3, bottom panel) by 95 Ma (Wynne et al., 1995; Enkin et al., 2003). Upper Cretaceous sedimentary and volcanic strata of the Silverquick–Powell Creek succession were then subsequently deposited unconformably on the combined Intermontane-Insular crustal block in this southerly position between 95 and 85 Ma, and the entire crustal block was subsequently translated ∼3000 km northward to its present position by ca. 52 Ma, a complicated journey at best (Fig. 3, bottom panel; Enkin, 2006).

The high latitudinal translation rates required by the paleomagnetic data are also problematic and contrast markedly with those interpreted from Pacific plate kinematics (Engebretson et al., 1985). The ∼2000 km of southward motion required for the Intermontane-Insular block after accumulation of the Spences Bridge volcanic rocks would have had to take place during a period (100–95 Ma) when the Farallon and Kula plates were thought to have been subducting orthogonally (Fig. 3, bottom panel; Stock and Molnar, 1988). The inferred translation rates are rapid for such a large crustal block, even if it were coupled directly to the oceanic plate (Gordon, 1998). In addition, translation of this continental-scale block must have taken place while leaving no structural or stratigraphic evidence of its rapid movement along, and interaction with, the continental margin, despite having had to pass by the Sierra Nevada magmatic arc twice in Late Cretaceous–Cenozoic time. There is no evidence for the development of synkinematic transtensional basins along the continental margin, nor evidence for marginal basins or ocean-floor material accreted along a Cenozoic suture on the eastern side of the combined Intermontane-Insular crustal block. In fact, there is no geological or geophysical evidence in the western North American Cordillera for any crustal-scale translational structure(s) that could have accommodated large-scale sinistral translation (∼2000 km) followed by dextral translation (∼3000 km) movements along the eastern edge of the combined superterranes that are required by the paleomagnetic data (Wyld et al., 2006; McMechan et al., 2021). The lack of these features must be accounted for in any Late Cretaceous–Cenozoic tectonic reconstruction of the Canadian Cordillera. Alternative models have been proposed, but none address the apparent lack of geologic evidence that would support large-scale translation (Johnston, 2008; Kent and Irving, 2010; Hildebrand, 2013).

Recent reconstructions (Matthews et al., 2017; Sauer et al., 2017a, 2018, 2019) focus on only limited portions of the paleomagnetic data and neglect to address the entirety of the data set; however, when considered in total, the full paleomagnetic data set requires that “Baja B.C.” must have been far larger than is commonly assumed and would have had to undergo multistage southward and northward latitudinal displacements, an extremely complicated tectonic history with no supporting geologic evidence.

The Nanaimo Basin provides a high-fidelity record of orogenic exhumation and basin evolution in the southwestern Canadian Cordillera that constrains the tectonic evolution of the combined Intermontane-Insular block. Analysis of the provenance and evolution of the Nanaimo Basin will therefore provide data crucial to resolution of the terrane translation issue.

Basin Nomenclature

Previous workers have referred to the entire Upper Cretaceous to Cenozoic tectonostratigraphic succession underlying and flanking the modern Georgia Strait as the Georgia Basin (England, 1989; Mustard, 1994). This nomenclature belies the complicated depositional history of this tectonostratigraphic succession. The succession consists of at least two unconformity-bounded stratigraphic sequences, including the largely marine Upper Cretaceous–Paleocene Nanaimo Group and the overlying nonmarine Paleocene to Oligocene Huntingdon and Chuckanut formations (Johnson, 1984; England and Bustin, 1998) and temporally equivalent marginal marine rocks (Haggart et al., 2018a). These unconformity-bounded sequences differ fundamentally in age, provenance, depositional environments, stratigraphic evolution, and structural setting and were not deposited in a single depocenter. It is more parsimonious to recognize that the Georgia Basin is a northwest-southeast–oriented neotectonic structural and topographic depression that contains two distinct stratigraphic successions, including strata of both the Late Cretaceous to Paleocene Nanaimo Basin (Pacht, 1984; Katnick and Mustard, 2003) and the unconformably overlying Paleocene to Oligocene Chuckanut Basin (Johnson, 1984).

Prior to recent studies, most workers assumed that deposition of the Nanaimo Group ended in latest Maastrichtian time, although some data suggested that Paleocene strata might be present locally in the succession (Muller and Jeletzky, 1970; Ward et al., 2012). The high percentage of syndepositional detrital zircon in much of the Nanaimo Group makes it an ideal candidate for the application of maximum depositional age (MDA) calculations in assessing the absolute depositional age of stratigraphic units. The recent analysis of large-n detrital zircon data sets provides high-resolution chronostratigraphic control, particularly within the upper portion of the Nanaimo Group (Englert et al., 2018, 2020; Coutts et al., 2020). MDAs from the non-fossiliferous Gabriola Formation in the northern Nanaimo Basin support the suggestion that sediment accumulation within the basin continued into the early Paleocene (ca. 63–64 Ma) (Matthews et al., 2017; Coutts et al., 2020). We follow these recent workers herein and consider that deposition of the Nanaimo Group succession continued into the earliest Paleocene.

The Nanaimo Group is exposed primarily along the eastern shore of Vancouver Island and in islands immediately offshore, with isolated exposures on the western coast of the British Columbia mainland (Fig. 2). The group consists of a >4 km thickness of clastic sedimentary strata of Turonian to Paleocene age that were deposited on the southeastern portion of the Wrangellia terrane of the Insular assemblage, on the western margin of an active continental magmatic arc currently represented by the dissected southern Coast Mountains batholith (Figs. 2A, 2B). The modern distribution of this group is the result of differential exhumation along both contractional and extensional structures (England, 1989; England and Bustin, 1998).

The precise basin classification of the Nanaimo Basin has been a matter of debate (Pacht, 1984; England, 1989; Mustard, 1994). The basin is flanked to the east by a Jura-Cretaceous magmatic arc and to the west by an outer-arc high represented by the Wrangellia terrane, and thus fits the classic definition of a forearc basin (Dickinson and Seely, 1979; Ward and Stanley, 1982; Dickinson, 1995). However, initial subsidence in the Nanaimo Basin was driven by west-vergent contractional deformation within the Coast Belt thrust system (Journeay and Friedman, 1993) and the North Cascades thrust system (Pacht, 1984; Brown, 2012), leading Mustard (1994) to characterize the Nanaimo Basin as a foreland basin. A modern classification scheme based on material transfer between the subducting and overriding plate and the style of deformation within the upper plate suggests that the Nanaimo Basin was initiated as a “compressional accretionary forearc basin” (Noda, 2016). This classification is supported by evidence of: (1) underplating of Nanaimo Basin strata under the North Cascades (Matzel et al., 2004; Gordon et al., 2017; Sauer et al., 2017b); (2) contractional structures in the adjacent arc massif to the east (Coast Belt and North Cascades thrust systems) (Journeay and Friedman, 1993; Brown, 2012); and (3) initial extensional(?) uplift of the outer-arc high in Turonian–Santonian time (ca. 93–85 Ma) (England and Bustin, 1998).

Biostratigraphic Age Constraints

The age range of marine strata within the Nanaimo Basin, representing arguably >90% of the succession, is well constrained biostratigraphically. Muller and Jeletzky (1970) established the modern biostratigraphic subdivisions of the Nanaimo Group sensu stricto based on ammonites and inoceramid bivalves, and this framework has been subsequently refined and improved (Ward, 1978; Haggart et al., 2005, 2009; Ward et al., 2012; Haggart and Graham, 2018). This modern biostratigraphic framework, supported by microfossil studies (McGugan, 1962), has demonstrated that the typical succession of the Nanaimo Group ranges in age from mid-Santonian to Maastrichtian. Stratigraphic relations at the base of the Nanaimo Group succession, however, are problematic. Turonian and Coniacian ammonites and inoceramid bivalves have been identified in strata mapped at the base of the succession (Haggart, 1991a, 1994; Haggart et al., 2003), and Haggart et al. (2005) argued that these Turonian and Coniacian strata are stratigraphically distinct from the mid-Santonian Comox Formation that forms the base of the Nanaimo Group in most areas (Muller and Jeletzky, 1970). Unraveling the stratigraphic complexities within the lower part of the Nanaimo Group succession is strongly dependent on biostratigraphic data and requires the use of faunas that are recognized to have global correlation value.

Recent investigations have utilized MDAs of detrital zircon to constrain the approximate depositional age of Nanaimo Group strata (Matthews et al., 2017; Huang et al., 2019; Coutts et al., 2020). This approach provides excellent age constraints in strata that are biostratigraphically well constrained. However, MDAs of samples from strata that have been mapped previously as Comox Formation but that lack biostratigraphic constraints have been used to suggest that significant disconformities exist within the lower Nanaimo Group and that the base of the group may be as old as Middle Jurassic or Early Cretaceous (Huang et al., 2019). The presence of Late Cretaceous zircon in a sample yielding a Middle Jurassic MDA indicates that the quoted MDAs do not represent true depositional ages (Huang et al., 2019; Copeland, 2020; Kent et al., 2020).

Lithostratigraphic Succession

The Nanaimo Group consists of multiple clastic units of conglomerate, sandstone, and mudstone, with appreciable coal in some formations (Figs. 4, 5). The majority of the group consists of marine sandstone-conglomerate successions separated by thick intervals of mudstone and fine- to medium-grained sandstone. Formation boundaries are conformable and gradational, and contacts commonly interfinger on a scale of meters to kilometers. Internal lateral and vertical thickness changes and facies variations are common, and abrupt facies migrations have produced locally sharp and erosive formation boundaries (Figs. 4, 5). Sandstone intervals of the Nanaimo Group are dominated by submature, moderately sorted, medium- to coarse-grained feldspathic arenite, but there are distinct lateral and vertical variations in petrography (Pacht, 1984; Mustard, 1994; Coutts et al., 2020). The lower portion of the Nanaimo Group contains a higher percentage of chert and volcanic lithic fragments and is petrographically heterogeneous compared to the upper portion of the Nanaimo Group. Chert lithic arenite is present only in the southern portion of the Nanaimo Basin, primarily in the lower Nanaimo Group (Pacht, 1984; Coutts et al., 2020). The lower Nanaimo Group throughout the basin also contains a higher percentage of volcanic lithic fragments, although these are also noted locally in the upper Nanaimo Group (Pacht, 1984; Coutts et al., 2020). The upper Nanaimo Group is a distinctly more homogeneous package of feldspathic arenite. A distinct increase in plagioclase and potassium feldspar and a decrease in lithic fragments has been documented within the upper Nanaimo Group (Pacht, 1984; Mustard, 1994). This increased feldspar component corresponds to the appearance of locally conspicuous medium- to coarse-grained detrital muscovite in the upper Nanaimo Group.

Differences in coal stratigraphy and micropaleontology in the lower Nanaimo Group as well as physical separation of the basin into two depocenters by a topographic high (Nanoose arch) led earlier workers to subdivide the outcrop areas into a southern Nanaimo subbasin and northern Comox subbasin (McGugan, 1962; Bickford and Kenyon, 1987). Subsequent research has suggested that these stratigraphic differences in part reflect variable sedimentation along a topographically complex unconformity at the base of the succession (Muller and Jeletzky, 1970; Haggart, 1991a, 1994; Mustard, 1994; Haggart et al., 2005; Johnstone et al., 2006; Jones et al., 2018; Kent et al., 2020). In addition, the validity of the subbasin designation is questionable because demonstrable lateral continuity of strata along the length of the basin by early Campanian time indicates that the Nanaimo Basin was a large contiguous depocenter, and there is thus little utility in maintaining separate basin nomenclature after that time (Muller and Jeletzky, 1970; Mustard, 1994; England and Bustin, 1998).

The basal unconformity of the Nanaimo Basin is highly irregular, with paleo-relief locally >100 m (Mustard, 1994; Johnstone et al., 2006). Turonian-, Coniacian-, and Santonian-age deposits of the succession unconformably overlie Devonian metavolcanic rocks, Triassic basalt, and Jurassic plutonic rocks of the Wrangellia terrane (Muller and Jeletzky, 1970; Haggart, 1991a; Mustard, 1994; Haggart et al., 2005). The oldest stratigraphic units recognized in the basin are of Turonian and Coniacian age and are geographically restricted; they have not all been formally defined. These deposits resulted from initial subsidence of the basin in early Turonian time in response to contractional crustal thickening in the Coast Mountains batholith and Cascade Range to the east/southeast (Haggart, 1991a; Journeay and Friedman,1993; Haggart et al., 2005; Brown, 2012). Haggart et al. (2005) argued that these older Turonian–Coniacian strata represent a separate cycle of sedimentation prior to deposition of the younger Comox Formation, and, as such, they comprise distinct stratigraphic successions exposed only locally at the southern and northern limits of the Nanaimo Basin (Figs. 4, 5).

The first stratigraphic unit in the Nanaimo Group succession with basin-wide extent is the Comox Formation, a heterogeneous assemblage of conglomerate, sandstone, siltstone, and shale with local coal measures, biostratigraphically dated as Santonian in age (Figs. 4, 5; Muller and Jeletzky, 1970). Early Campanian basal transgressive equivalents of the Comox Formation have been mapped locally on Nanoose arch (Haggart et al., 2011), indicating the diachronous nature of initial marine sedimentation as it onlapped topographic highs in the basin. A diverse array of depositional environments is reflected in the Comox Formation, from nonmarine alluvial fan, braided river, and coastal plain environments to nearshore shallow-marine facies reflecting storm-dominated rocky shorelines and upper- to lower-shoreface conditions (Mustard, 1994; Haggart et al., 2005; Johnstone et al., 2006; Jones et al., 2018; Kent et al., 2020). Paleocurrent and lithofacies analyses suggest sediments of the Comox Formation were derived primarily from the south and southwest (England, 1989; Mustard, 1994).

The Comox Formation is conformably overlain by mudstone of the mid- to outer-shelf Haslam Formation, itself overlain by marine sandstone and conglomerate of the Extension and Protection formations. This succession of stratigraphic units, comprising the Comox through Extension and Protection formations, as well as older unnamed strata of Turonian and Coniacian ages, constitutes an overall transgressive sequence that progressed from nonmarine and coastal-plain strata through marginal-marine and marine strata of the Haslam and Protection formations, and is referred to informally as the lower Nanaimo Group succession (Figs. 4, 5, 6; England, 1989; Haggart, 1991a; Mustard, 1994; Kent et al., 2020). Jones et al. (2018) attributed the development of a transgressive, deepening-upward succession at the base of the Comox Formation in early Santonian time in the southern Nanaimo Basin to tectonic subsidence and subsequent rapid creation of accommodation space. Kent et al. (2020) applied a sequence stratigraphic framework to the lower succession in the northern Nanaimo Basin, documenting seven distinct depositional phases within an overall transgressive sequence.

In contrast to the overall nonmarine, marginal-marine, and shallow-marine nature of the lower Nanaimo Group succession, strata of the upper succession of the Nanaimo Group comprise deep-water (sub–wave base), fully marine lithofacies of early Campanian to Paleocene age (Figs. 4, 5, 6). The upper succession of the Nanaimo Group is composed of >2.1 km of deep-marine strata (Mustard, 1994; Katnick and Mustard, 2003) and represents more than two-thirds of the thickness of the group. The succession consists of six formations that alternate between recessive, laterally continuous, thin-bedded siltstone- and mudstone-dominated intervals (Cedar District, Northumberland, Spray formations) and thick-bedded to massive, lenticular conglomerate and sandstone intervals (De Courcy, Geoffrey, Gabriola formations) (Figs. 4, 5, 6). The formations are stacked vertically and interfinger laterally on a scale of tens to thousands of meters. The basal contacts of the coarse-grained units tend to be sharp and erosive, and the upper contacts gradational into overlying fine-grained units (Katnick and Mustard, 2003; Bain and Hubbard, 2016). Contacts between formations are regionally conformable and represent lateral and vertical migration of depositional systems, not major unconformable surfaces (Katnick and Mustard, 2003).

Bain and Hubbard (2016) conducted a detailed facies analysis of the upper succession of the Nanaimo Group on Denman and Hornby islands (Fig. 2B), subdivided the strata into ten lithofacies, and outlined a stratigraphic architecture that demonstrates convincingly that these strata represent a long-lived (∼15 m.y.) composite submarine channel system (20 km wide, 1.5 km thick) forming a slope conduit into a major submarine fan complex. Paleoflow measurements document southwest-directed paleoflow (227°; n = 2089), consistent with sediment derivation from the Coast Mountains batholith to the east (Bain and Hubbard, 2016). The lateral and vertical extent of this channel system indicates that the upper succession of the Nanaimo Group originally extended much further to the west and that the current outcrop belt is the erosional remnant of an initially much larger system (England, 1989: Haggart, 1993; Mustard, 1994; Kent et al., 2020; Englert et al., 2018).

The transition between the lower and upper successions of the Nanaimo Group took place rapidly across the basin during early Campanian time (Ward and Stanley, 1982; Mustard, 1994; Haggart et al., 2005; Katnick and Mustard, 2003; Bain and Hubbard, 2016), leading to progradation of submarine-fan systems to the west/northwest across the basin and subsequent accumulation of a thick complex (>2 km) of intertonguing fan lobes. Lateral and vertical variations in stratigraphic architecture indicate a multiphase submarine-channel system evolution that progressed from early sediment bypass and incision, to lateral channel migration, to a period of sediment aggradation, deep incision (>500 m), and bypass, and a final phase of aggradation (Bain and Hubbard, 2016; Englert et al., 2018, 2020). Detailed facies analysis combined with detrital zircon geochronology has permitted high-resolution stratigraphic reconstruction of the upper succession of the Nanaimo Group, providing excellent constraints on deep-water sediment dispersal within a Late Cretaceous forearc basin system (Bain and Hubbard, 2016; Englert et al., 2018, 2020).

The upper succession of the Nanaimo Group is overlain by an early Paleocene to Eocene succession of marginal marine to nonmarine conglomerate, sandstone, and siltstone (Johnson, 1984; England and Bustin, 1998; Mustoe et al., 2007; Matthews et al., 2017). These strata, including the Huntingdon and Chuckanut formations, are exposed primarily in the southern and southeastern portion of the Georgia Basin depocenter. Scattered outcrops of presently undescribed strata of Paleocene age are preserved on the eastern side of Vancouver Island and on the southwestern Gulf Islands (Fig. 2; Haggart et al., 2018a). Uplift and erosion of the Nanaimo Group in post–early Paleocene time is suggested by a dramatic change in depositional environment from deep water to marginal marine, locally significant relief at the base of the overlying succession, and vitrinite reflectance data that suggest removal of as much as 1000 m of probable Maastrichtian- to Paleocene-age strata (England and Bustin, 1998).

High-resolution provenance analysis is a prerequisite for accurate reconstruction of the tectonic evolution and paleogeographic setting of the Nanaimo Basin. Overlap among detrital zircon age populations from multiple proposed source regions leads to nonunique solutions (Matthews et al., 2017), which exemplifies the importance of multidisciplinary provenance analysis. This investigation integrates detrital zircon U-Pb geochronology, conglomerate clast U-Pb geochronology, detrital muscovite 40Ar/39Ar thermochronology, and zircon Lu-Hf isotopic analysis of detrital zircon to constrain the source region(s) supplying sediment to the Nanaimo Basin in Late Cretaceous time.

Detrital Zircon U-Pb Geochronology

Multiple investigations have examined detrital zircon U-Pb geochronology of the Nanaimo Group since the mid-1990s. Interestingly, these studies precisely track the dramatic evolution in analytical methodology in detrital zircon geochronology over the past 25 years (e.g., Gehrels et al., 1995; Gehrels, 2011, 2014; Gehrels and Pecha, 2014). Mustard et al. (1995) utilized single-crystal U-Pb thermal ionization mass spectrometry (TIMS) to examine detrital zircon selected by morphology and color and were the first to recognize sediment derivation from multiple sources, including Proterozoic metasedimentary rocks and Mesozoic plutons. Sensitive high-resolution ion microprobe (SHRIMP) analysis of bulk zircon separates (n ≈ 60) expanded the number and stratigraphic extent of detrital zircon ages in the Nanaimo Group and verified the presence of both local and extra-regional sediment sources (Mahoney et al., 1999, 2003; Mustard et al., 2007). Laser ablation inductively coupled plasma mass spectrometry (LA-ICPMS) increased both the number of zircons per sample (n = 60–200) and the geographic extent of detrital zircon analyses and confirmed the regional consistency of detrital zircon patterns (Mahoney et al., 2014, 2016; this study). Recently, large-n (>300) detrital zircon geochronology has been utilized to constrain the depositional age and provenance of the Nanaimo Group (Matthews et al., 2017; Huang, 2018; Huang et al., 2019; Coutts et al., 2020; Kent et al., 2020). These analyses have substantially increased the precision of the detrital zircon data set and have confirmed detrital zircon patterns established by earlier investigations. Significantly, these large-n detrital zircon analyses have yielded high-precision determination of MDAs that permit detailed reconstructions of stratigraphic architecture (Haggart et al., 2018b; Huang et al., 2019; Englert et al., 2018, 2020; Coutts et al., 2020; Kent et al., 2020).

This investigation reports the results of detrital zircon analyses (25 samples) from three composite transects across the north, central, and southern portions the Nanaimo Basin (Figs. 2B, 5, 79; Items A and B in the Supplemental Material1). The northern transect encompasses Vancouver, Hornby, and Denman islands and is the most stratigraphically complete sample suite (14 samples from eight formations) (Figs. 2B, 4, 5, 7). The central transect includes primarily the Nanaimo area and Gabriola Island and is more limited in stratigraphic extent (four samples from four formations). This second transect includes legacy TIMS data and two samples containing Archean zircons previously described in Mahoney et al. (1999). The southern transect encompasses Galiano and Salt Spring islands and focuses on the coarse-grained lithofacies in the section (seven samples from six formations) (Figs. 2B, 4, 5, 7). Two samples from the eastern basin margin, one from strata equivalent to the Nanaimo Group and one from Eocene strata of the Huntingdon Formation, were also analyzed. Detrital zircon samples were collected from localities with associated biostratigraphic control to constrain stratigraphic age (see Supplemental Items A, B). Detrital zircon data are presented as probability density plots in Figure 7, and analytical data are presented in Supplemental Item C. Multiple analytical methods (TIMS, SHRIMP, and LA-ICPMS) were employed over the course of this investigation; analytical details are provided in Supplemental Item C.

Conglomerate Clast U-Pb Geochronology

Conglomerates are excellent provenance indicators because they are deposited proximal to their source, recording hundreds of kilometers or less, not thousands of kilometers, of transport (Paola, 1988). Quartzite clasts are particularly useful because the high concentration of zircon in these clasts allows U-Pb geochronology to provide a high degree of confidence about clast provenance.

Conglomeratic facies are a conspicuous component of the coarse-grained intervals of the upper Nanaimo Group and are particularly well developed in the Geoffrey Formation (Mustard, 1994; Katnick and Mustard, 2003; Bain and Hubbard, 2016; Englert et al., 2018, 2020). Conglomerate clast lithologies are dominated by volcanic rocks, intermediate plutonic rocks, and chert, with lesser mafic volcanic, metamorphic, and sedimentary rocks (Mustard, 1994; Katnick and Mustard, 2003). A small yet conspicuous component of the conglomerate lithofacies consists of white and pink, pebble- to boulder-sized quartzite clasts (Fig. 10).

Five quartzite clasts were collected from the Geoffrey Formation, including four on Hornby Island and one on Salt Spring Island (Figs. 5, 10, 11). Detrital zircon data are presented as probability density plots in Figure 11, and analytical data are presented in Supplemental Item D. Three of the clasts were analyzed by SHRIMP and two of the clasts were analyzed by LA-ICPMS; analytical details are provided in Supplemental Item D.

Detrital Muscovite 40Ar/39Ar Data

Detrital muscovite forms a minor, yet conspicuous, component of the upper Nanaimo Group, particularly in the Gabriola Formation, where millimeter-sized flakes may be found on bedding planes. Detrital muscovite was separated from one sample of the middle Campanian De Courcy Formation on Salt Spring Island and from two samples of the Maastrichtian Gabriola Formation (lower and upper portions) on Galiano Island (Figs. 5, 12). These samples were processed for 40Ar/39Ar thermochronology at the Pacific Centre for Isotopic and Geochemical Research at the University of British Columbia (Vancouver, British Columbia, Canada). Analytical methodology and isotopic data are described in Supplemental Item E.

Detrital Zircon Hafnium (Hf) Isotopic Geochemistry

Combined detrital zircon U-Pb geochronology and Lu-Hf analysis were conducted on five samples from the northern Nanaimo Basin. Subpopulations of previously dated grains were selected from the Cedar District (163–93 Ma), De Courcy (96–74 Ma), Geoffrey (81–73 Ma), Spray (106–65 Ma), and Gabriola (78–66 Ma) formations in order to characterize Mesozoic age subpopulations through the upper Nanaimo Group (Figs. 5, 13). The Lu-Hf isotopic values were obtained from the LaserChron Laboratory at the University of Arizona (Tucson, Arizona, USA). Grains were selected to represent each of the main age subpopulations and to avoid crystals with discordant or imprecise ages. Laser spots were placed over the U-Pb analysis pits in order to ensure that Hf isotopic data were determined from the same domain as the U-Pb age. Hf data are presented on an Hf-evolution diagram (Fig. 13; Supplemental Item F).

Provenance Constraints: Detrital Zircon Data

The bimodal character of the detrital zircon distribution within the Nanaimo Group, consisting primarily of Mesozoic (ca. 160–70 Ma) and Proterozoic (ca. 1800–1350 Ma) populations, has long been recognized (Mustard, 1994; Mustard et al., 2007; Mahoney et al., 1999, 2014, 2016; Matthews et al., 2017; Coutts et al., 2020). In this study, the Mesozoic and Proterozoic detrital zircon populations are discussed separately in order to develop a high-fidelity assessment of provenance variations with the Nanaimo Group (Figs. 7, 8, 9).

Mesozoic Detrital Zircon Population

Approximately 81% (n= 1341/1663) of detrital zircons analyzed in this investigation are of Phanerozoic age, and of these, 94% (n= 1263/1341) are Mesozoic in age (Figs. 7, 8; Supplemental Item C). Probability density plots, normalized age probability plots, and cumulative distribution functions (Gehrels, 2011; Saylor and Sundell, 2016) were utilized to analyze the lateral and vertical detrital zircon distribution within the Nanaimo Group. These analyses permit subdivision of the Mesozoic zircon population into four distinct detrital zircon populations (DZPs) (Fig. 8B).

Subpopulation DZP1 is characterized by a distinct Devonian–Mississippian population, with a ca. 353 Ma (380–330 Ma) peak, and a broad trimodal Jurassic–Cretaceous population, with peaks at ca. 97 Ma (105–85 Ma), ca. 123 Ma (130–115 Ma), and ca. 155 Ma (175–140 Ma) (Fig. 8). The DZP1 signature is restricted to unnamed strata of early Turonian to Coniacian age in the northern part of the basin and the late Santonian to early Campanian Extension Formation in the northern and southern transects (Fig. 5). These formations were not sampled in the central transect. Proterozoic grains are virtually absent (0.6%) from DZP1. The Devonian–Mississippian age population is restricted to the southern Nanaimo Basin and is absent in the north.

Subpopulation DZP2 is characterized by a bimodal Jurassic–Cretaceous population, with peaks at ca. 150 Ma (165–140 Ma) and ca. 88 Ma (100–80 Ma), and a subordinate ca. 75 Ma (80–73 Ma) peak (Fig. 8). The DZP2 signature derives from the early to middle Campanian Protection (two samples), Cedar District (three samples), and De Courcy (four samples) formations. Proterozoic grains comprise 22% of DZP2 grains but are present only in the southern and central transects.

Subpopulation DZP3 has a very distinct unimodal peak at ca. 73 Ma (82–69 Ma) that encompasses 56% of all grains in the sample suite (Fig. 8). The DZP3 signature is aggregated from the upper Campanian Northumberland (one sample), uppermost Campanian to Maastrichtian Geoffrey (two samples), and Maastrichtian–Paleocene Gabriola (one sample) formations (Figs. 7, 8). Proterozoic grains represent 19% of DZP3 zircon grains and are found in all three transects.

Subpopulation DZP4 displays a broad, trimodal Jurassic–Cretaceous population with peaks at ca. 70 Ma (75–65 Ma), ca. 89 Ma (105–80 Ma), and ca. 161 Ma (190–145 Ma) (Fig. 8). The DZP4 signature is derived from the Maastrichtian Geoffrey (n = 2) and Maastrichtian–Paleocene Spray (n = 1) and Gabriola (n = 4) formations from the northern and central transects (Fig. 7). Proterozoic grains represent 29% of all DZP4 zircon grains. The Gabriola Formation is the stratigraphically highest formation in the Nanaimo Group. The youngest MDA, calculated by the YGC2σ (youngest grain cluster at 2σ uncertainty) method of Dickinson and Gehrels (2009), from the Gabriola Formation in this study is 67.9 ± 1.6 Ma. Coutts et al. (2020) reported an MDA of 63.9 ± 1.9 Ma from correlative strata on Hornby Island, which is assumed to be more accurate due to the large-n analytical method.

The DZPs defined here are lithostratigraphically specific and characteristic of different formations within the Nanaimo Group. The formations interfinger laterally across the basin and are therefore commonly diachronous, which means the DZPs are also diachronous. For example, in the upper portion of the upper Nanaimo Group, DZP3 is found in the Northumberland and Geoffrey formations in the northern transect and in the Gabriola Formation in the central transect. This suggests that DZP3 is diachronous from north/northeast to south/southeast, becoming progressively younger to the south/southeast.

Proterozoic Detrital Zircon Population

Approximately 19% of the detrital zircons analyzed in this investigation are >500 Ma. The vast majority are Proterozoic (98%) in age. Archean grains compose ∼2% of the data set, although Coutts et al. (2020) reported a substantially lower percentage (0.2%) using a large-n analytical approach. Less than 1% of the grains are Cambrian in age. The Proterozoic component of the data set was subdivided into the same DZP bins defined in the Mesozoic component to capture lateral and vertical variations in provenance within the Nanaimo Basin (Fig. 9). Data are displayed as an aggregate probability density plot (Fig. 9A) and as a series of normalized age probability plots (Figs. 7, 9B). DZP1 contains two zircon grains of Proterozoic age (0.06%; n = 2/338) and is not included in Figure 9.

Subpopulation DZP2 displays a broad bimodal zircon distribution with peaks at 1368 Ma (1415–1335 Ma) and 1702 Ma (1800–1650), with minor early Proterozoic–late Archean (2670–2240 Ma) and Grenville-age (1100–900 Ma) populations (Fig. 9). There is also a small bimodal peak between 625 Ma and 585 Ma. Importantly, the Proterozoic signal is not synchronous across the basin but is restricted to the late Santonian to lower Campanian Protection and middle Campanian De Courcy formations in the southern and central transects (Fig. 7). Proterozoic zircon is not evident in the Cedar District and De Courcy formations in the northern transect and did not reach the northern portion of the basin until Maastrichtian time, clearly illustrating the diachronous nature of sedimentation in the basin.

Subpopulations DZP3 and DZP4 are statistically indistinguishable (Kolmogorov-Smirnov [K-S] p-value = 0.528) with sharp bimodal peaks at 1707 and 1719 Ma, respectively (1790–1690 Ma), and 1386 and 1376 Ma, respectively (1430–1325 Ma). Both the Paleoproterozoic–late Archean and Grenville-aged populations are significantly reduced. There is a scattering of 700–500 Ma grains in both subpopulations. The bimodal detrital zircon distribution is characteristic of the upper two-thirds of the upper Nanaimo Group (Northumberland, Geoffrey, Spray, and Gabriola formations) across the basin. The only exception to the ubiquitous nature of this detrital zircon signature in the upper Nanaimo Group is in the Gabriola Formation within the central transect, which displays a significantly diminished Proterozoic signature (Fig. 7).

Conglomerate Clast Data

Each of the five quartzite clasts collected from the Geoffrey Formation has a unique detrital zircon U-Pb spectrum that is statistically distinct from that of the other clasts (Fig. 11). There are broad similarities between the quartzite clast detrital zircon spectra, including the presence of a strong unimodal 1800–1700 Ma peak and a small number of Paleoproterozoic to Archean grains. Differences include the age range of a minor 1450–1300 Ma peak, the presence or absence of Grenville-aged zircon, and the presence or absence of non–North American magmatic gap (1610–1490 Ma) grains (Fig. 11; Supplemental Item D).

Four of the five quartzite clasts collected from the Geoffrey Formation have zircon U-Pb age distributions that statistically match with potential source regions in British Columbia, Montana, and Idaho based on the statistical comparisons (i.e., K-S test; Figs. 10, 11). Interestingly, each of the clasts appears to be derived from different Mesoproterozoic to Paleogene geologic units. These include: (1) the Phillips Formation of the Mesoproterozoic Purcell Supergroup in southern British Columbia (Gardner, 2008); (2) the Apple Creek Formation of the Mesoproterozoic Lemhi subbasin in central Idaho (Burmester et al., 2016); (3) the Cambrian Flathead Sandstone of west-central Montana; and (4) quartzite clasts from the Paleocene Pinyon Formation in southwestern Montana, which are interpreted to be derived from the Neoproterozoic Brigham Group (Malone et al., 2017). These multiple high-resolution correlations link the conglomeratic strata of the Geoffrey Formation with Mesoproterozoic to Cambrian strata within and overlying the Belt-Purcell Basin of Montana, Idaho, and British Columbia (Figs. 10, 11).

Detrital Muscovite 40Ar/39Ar Data

Detrital muscovite age spectra are displayed in Figure 12. Muscovite from the base of the upper Campanian De Courcy Formation (sample 14JBM06) in the southern Nanaimo Basin yielded an 40Ar/39Ar age of 78.33 ± 0.48 Ma (Fig. 12). Muscovite from the basal portion of the Maastrichtian Gabriola Formation in the southern Nanaimo Basin (sample 16JBM06) yielded an 40Ar/39Ar age of 76.70 ± 0.42 Ma, and a sample from the upper portion of the formation (sample 17JBM06) yielded an overlapping 40Ar/39Ar age of 77.32 ± 0.44 Ma (Fig. 12). These results demonstrate that the increase of young (ca. 79–73 Ma) zircon in DZP2 was accompanied by an influx of essentially coeval muscovite, representing a distinct pulse of syndepositional detritus into the upper Nanaimo Basin.

Detrital Zircon Hafnium (Hf) Isotope Data

Detrital zircon Hf data from the upper Nanaimo Group fall into three distinct clusters (Fig. 13). The first cluster is characterized by positive εHf values (>+9) and contains all zircon older than 85 Ma, including all samples from the Cedar District and De Courcy formations. This cluster also includes a small subset of young zircon (75–65 Ma) from the De Courcy and Spray formations. The second cluster displays “transitional” εHf values (∼+5 to −2) and characterizes subpopulations of young zircon (ca. 80–65 Ma) from the Geoffrey, Spray, and De Courcy formations (Fig. 13). The third cluster is characterized by negative εHf values (>−10) and includes young zircon (78–65 Ma) from the Geoffrey, Spray, and Gabriola formations (Fig. 13).

Detrital Zircon Trends

The DZPs defined in the Nanaimo Group constitute discrete provenance signatures that constrain the depositional history of the basin. Subsidence in the basin began in the early Turonian, coincident with contraction in the southern Coast Mountains batholith that uplifted Middle Jurassic to Early Cretaceous plutons and pre- and synkinematic plutons to provide the broad Jurassic–Cretaceous trimodal population evident in DZP1 (Figs. 7, 8). The ca. 123 Ma peak is distinctive in that there is a well-defined magmatic gap in the batholith from 120 to 140 Ma (Gehrels et al., 2009). However, zircon of this age has been documented in the Methow Basin, and its occurrence in the lower Nanaimo Group is consistent with inversion of that basin during regional contraction in Cenomanian to Turonian time (DeGraaf-Surpless et al., 2003; Mahoney et al., 2014). The distinct Devonian–Mississippian component in DZP1 is restricted to the southern Nanaimo Basin and is assumed to be locally derived from the Late Devonian Salt Spring Island intrusions and the mid-Paleozoic Sicker Group (Figs. 7, 8).

There was a significant shift in provenance in late Santonian time (ca. 85–84 Ma), coincident with a major increase in subsidence, beginning with the deposition of the Protection Formation in the southern Nanaimo Basin. The provenance shift is characterized by a transition to a younger Mesozoic signature, with a bimodal peak at ca. 88 Ma and ca. 75 Ma, as well as a major influx of Proterozoic zircon (Figs. 7, 8, 9). The ca. 75 Ma peak is weak in the south and becomes much more prominent upsection to the north in middle Campanian strata of the Cedar District and De Courcy formations. The diachronous nature of the provenance shift is particularly evident in the Proterozoic signal, which is prominent in the upper Santonian Protection Formation in the south but did not arrive in the north until the upper Campanian Northumberland Formation (Fig. 7).

There was a second substantial change in provenance and basin dynamics in the late Campanian–early Maastrichtian (ca. 71–70 Ma) involving a doubling of sedimentation rate accompanied by a flood of syndepositional detrital zircon, with 51% of DZP3 grains between 75 and 68 Ma (Figs. 7, 9, 14, 15). The very short lag time between crystallization and deposition of the 75–68 Ma zircon requires rapid exhumation of the source plutons. This flood of syndepositional zircon is first evident in the northern Nanaimo Basin during deposition of the Northumberland Formation and continues through the deposition of the Geoffrey Formation accompanied by a major pulse of channelized pebble to boulder conglomerate (Mustard, 1994; Englert et al., 2018, 2020; Coutts et al., 2020). The influx of syndepositional zircon is synchronous with a pulse of Proterozoic zircon into the northern Nanaimo Basin (Figs. 7, 9; Coutts et al., 2020). The Proterozoic zircon signal is strikingly uniform across the basin in Maastrichtian to Paleocene time, characterized by a distinct bimodal age spectrum with peaks at ca. 1704 Ma and ca. 1380 Ma.

There is a conspicuous increase in coarse-grained muscovite in the Gabriola Formation, accompanied by coeval (79–65 Ma) zircon, which requires rapid, syndepositional exhumation of Late Cretaceous muscovite-bearing plutons in the source region. The influx of Proterozoic zircon continued unabated throughout deposition of the Gabriola Formation.

A third change in provenance and basin dynamics is evident in the late Maastrichtian to early Paleocene (ca. 66–64 Ma), with an increase in both Late Cretaceous (ca. 89 Ma) and Jurassic (ca. 161 Ma) zircon (Fig 7, 8). This shift signals a reduction in the prominence of the syndepositional plutons and a return of older zircon from the southern Coast Mountains batholith within the Spray and Gabriola formations.

Detrital Zircon Provenance

The origin of the prominent ca. 75–68 Ma syndepositional detrital zircon population in the upper Nanaimo Group is not readily apparent. Late Cretaceous plutons are a minor component of the southeastern Coast Mountains batholith east of the Nanaimo Basin (Friedman and Armstrong, 1995; Cecil et al., 2018). Cecil et al. (2018) documented a voluminous magmatic event in the Knight and Bute Inlets to the northeast of the Nanaimo Basin at 85–70 Ma, which could have provided a source for much of the young zircon in the upper Nanaimo Group. Mahoney et al. (2009) recorded a major Late Cretaceous magmatic flareup in the Bella Coola region that could represent a northerly source. Late Cretaceous plutons are also abundant in the Idaho batholith, which could provide an extra-regional source to the east. Differentiation of potential sources of the ubiquitous young zircon in the upper Nanaimo Group requires Hf isotopic analyses, discussed below.

Proterozoic zircon in the upper Nanaimo Group is extra-regional (Mustard et al., 1995; Mahoney et al., 2003; Matthews et al., 2017). Proterozoic zircon is rare in Wrangellia, the Coast Mountains batholith, the Intermontane Belt, and the North Cascades. However, the Mesoproterozoic Belt-Purcell Supergroup of southern British Columbia, Idaho, and Montana is a prolific source of detrital Proterozoic zircon grains. The supergroup is composed of several different units, including the lower Belt, Ravalli, Pigeon, Missoula, and Lemhi groups, each of which displays a characteristic detrital zircon signature (Ross and Villeneuve, 2003; Link et al., 2006, 2016; Stewart et al., 2010). The Proterozoic age spectra in the Nanaimo Group are essentially bimodal, with strong peaks at ca. 1704 Ma and ca. 1380 Ma and minor associated populations of Archean, Paleoproterozoic, and Neoproterozoic grains (Fig. 9A). The narrow ca. 1704 Ma peak and scattered older zircon are very similar to those seen in the Lemhi subbasin of the upper Belt-Purcell Supergroup (Link et al., 2016). The Neoproterozoic population (ca. 600–500 Ma) may represent Ediacaran to Cryogenian plutons described from central Idaho (Lund et al., 2010). The lack of Grenville-age zircon in the Nanaimo Group is significant in that Neoproterozoic strata along the continental margin contain a substantial component of Grenville-age zircon (ca. 1100–900 Ma), suggesting that these strata were not a major sediment source for the Nanaimo Basin (Yonkee et al., 2014; Matthews et al., 2017).

It is important to note that 31% of the aggregate Proterozoic zircon population in the Nanaimo Basin is ca. 1380 Ma (1400–1320 Ma) in age (Fig. 9). Available detrital zircon data from Neoproterozoic to Eocene strata from Alaska to Mexico indicate that ca. 1380 Ma (1400–1320 Ma) zircon is extremely rare along the western Laurentian margin, accounting for <2% of all grains examined (Laskowski et al., 2013 [94 samples; n = 147/8717]; Gehrels and Pecha, 2014 [34 samples; n = 109/6164]). Zircon of this age is also rare in high-resolution studies of Neoproterozoic to Cambrian strata from central British Columbia to Mexico, representing <3% of all grains (Yonkee et al., 2014 [96 samples; n = 216/6785]; Matthews et al., 2017 [46 samples; n = 122/4780]). These regional data suggest that wherever the prolific source of ca. 1380 Ma zircon to the Nanaimo Basin was located, it was localized and not exposed prior to the Late Cretaceous and, once exposed, limited in the geographic dispersal of its zircon.

Dumitru et al. (2016) posited that the existence of a distinct ca. 1380 Ma age peak, particularly the presence of a “Lemhi doublet” (peaks at ca. 1380 Ma and 1800–1650 Ma), uniquely fingerprints strata of the Lemhi subbasin of the Belt-Purcell Supergroup and post-Belt plutonic rocks in central Idaho. The same unique signature has been identified in the northern Upper Cretaceous Great Valley forearc and Franciscan mélange, in the Eocene Tyee Group, and in the Upper Cretaceous Swakane Gneiss of the North Cascades (Dumitru et al., 2016). Lewis et al. (2010) reported ca. 1380 Ma detrital zircon in Cretaceous strata overlying ca. 1380 Ma gneiss in northern Idaho. Several ca. 1380 Ma plutons have been dated in central Idaho (Aleinikoff et al., 2012), and one of the authors of the current study (PKL) has analyzed detrital zircon from modern sediments of the Salmon River in central Idaho and documented a strong ca. 1380 Ma peak (∼20%, n = 11/52; Supplemental Item G). The Nanaimo Group contains the largest component of ca. 1380 Ma zircon found anywhere along the western continental margin, suggesting it was located proximal to its source. It is highly probable that the ca. 1380 Ma zircon in the Nanaimo Basin was derived from plutonic rocks in central Idaho.

Importantly, derivation of Nanaimo Group sediment from central Idaho is essentially required by detrital zircon U-Pb geochronology of quartzite clasts collected from boulder conglomerate of the Geoffrey Formation (Figs. 10, 11). These data demonstrate that these clasts were derived from Mesoproterozoic to Cambrian strata within and overlying the Belt Supergroup of Montana, Idaho, and British Columbia.

Detrital Muscovite 40Ar/39Ar Data

Muscovite is relatively rare in source regions immediately adjacent to the Nanaimo Basin. The Coast Mountains batholith is primarily tonalitic in composition, and muscovite is therefore sparse. Muscovite-bearing metamorphic rocks are currently exposed in the North Cascades, but these rocks were deeply buried in the Late Cretaceous and not exhumed until the Eocene (Matzel et al., 2004; Miller et al., 2016). Muscovite is locally abundant in the southeastern Coast Mountains batholith, west of the Fraser–Straight Creek fault, where it is found in metamorphic aureoles of Late Cretaceous plutons (95–84 Ma) and yields 40Ar/39Ar cooling ages of ca. 80 Ma (Brown et al., 2000; Mitrovic, 2013). However, these rocks were probably not exhumed until the Eocene, and this area has a paucity of plutons <84 Ma (Friedman and Armstrong, 1995; Cecil et al., 2018). While Paleocene to Eocene muscovite is common in the Omineca belt of southern British Columbia, in the Bitterroot lobe of the Idaho batholith, and in areas west of the Boulder batholith, Cretaceous muscovite is rare (Foster et al., 2001; Webster et al., 2020).

The Atlanta lobe of the Idaho batholith contains an extensive belt of two-mica granites produced between 83 and 68 Ma by substantial crustal thickening, crustal anatexis, and development of a large orogenic plateau (Fayon et al., 2017; Byerly et al., 2017). Slow exhumation of this orogenic plateau is consistent with central Idaho being a major sediment source for Late Cretaceous–Paleocene rivers in Washington, Montana, and California (Dumitru et al., 2016; Byerly et al., 2017). Fayon et al. (2017) analyzed a biotite granodiorite from the Atlanta lobe and utilized geochronology (U-Pb zircon 71.9 ± 2.7 Ma [T = 700–750 °C]) and thermochronology (40Ar/39Ar biotite, 69.0 ± 0.2 Ma [T = 300–400 °C]; [U-Th]/He zircon, 39.7 ± 5.3 Ma [T = 150 °C]) to model a two-stage cooling history for the batholith. The cooling history was characterized by initial rapid cooling and exhumation through the biotite closure temperature between 72 and 65 Ma, followed by slow cooling from 65 Ma to present. Thermochronologic constraints and the widespread occurrence of sediment derived from the Idaho batholith in Late Cretaceous–Eocene basins in California, Wyoming, Montana, and Oregon require that the batholith was topographically high and subject to >4 km of exhumation in Late Cretaceous to Paleocene time (Fayon et al., 2017). This degree of exhumation is consistent with the suggestion that the Atlanta lobe peraluminous suite may have originally had a larger footprint than it does today (Gaschnig et al., 2017). The existence of 83–68 Ma two-mica granites rich in both muscovite and zircon, coupled with evidence for a large orogenic plateau undergoing >4 km exhumation in Late Cretaceous to Paleocene time, makes the Atlanta lobe of the Idaho batholith an ideal source for voluminous sediment influx into Campanian to Paleocene strata of the Nanaimo Basin. K-feldspar 40Ar/39Ar thermochronology on samples from the northern Nanaimo Group transect and various basement sources in the Pacific Northwest by Isava et al. (2021) support this interpretation.

Zircon Hafnium (Hf) Data

Integration of detrital zircon age spectra and Lu-Hf isotopic data constrains potential source areas for strata of the Nanaimo Basin. Mesozoic zircon older than ca. 85 Ma is characterized by a restricted range of juvenile Hf signatures (εHf ∼+9 to +15) that are just below the depleted mantle value, indicating derivation from uncontaminated, mantle-derived crustal reservoirs. The age and isotopic character of these zircons are an excellent match for the uniformly primitive signature of plutonic rocks in the southern Coast Mountains batholith (Fig. 13; Friedman et al., 1995; Cecil et al., 2011).

The uniform, juvenile character of all zircon older than ca. 85 Ma is in marked contrast to the much broader range of Hf values in young zircon (80–65 Ma), which range from εHf ∼+10 to εHf ∼−30 (Fig. 13). Within this broad range, there are three distinct subpopulations, including juvenile (εHf = +10 to +15), transitional (εHf = −5 to +5), and evolved (εHf <−10), indicating a mixing of fundamentally different sources.

The juvenile values are consistent with derivation from Late Cretaceous plutons with uniformly juvenile Hf values that have been documented from the Knight Inlet–Bute Inlet region in the northern portion of the southern Coast Mountains batholith (Fig. 13A; Homan, 2017; Cecil et al., 2018).

The source of zircon with transitional Hf values is problematic. Hf and Nd isotopic data demonstrate that the vast majority of the Coast Mountains batholith is isotopically juvenile (CCMB and SCMB fields on Fig. 13A; Friedman et al., 1995; Cecil et al., 2011; Homan, 2017). The lowest Hf values (εHf ∼+2 to +5) are found on the western side of the northern batholith, where they are associated with mid-Cretaceous plutons intruding continental metasedimentary rocks of the Banks Island assemblage (Gehrels et al., 2009; Cecil et al., 2018). Plutons of Late Cretaceous age are not found along the western side of the batholith. It is plausible, but as yet unproven, that Late Cretaceous plutons intruding the continental-affinity Yukon-Tanana terrane within and along the eastern margin of the northern Coast Mountains batholith could yield transitional Hf signatures (Gehrels et al., 2009; Cecil et al., 2011; Pecha et al., 2016). There is a belt of Late Cretaceous plutons along the eastern margin of the southern Coast Mountains batholith that yield transitional Hf values (εHf ∼+4 to +5), but more isotopically evolved plutons have not been recognized in the southern batholith.

Surpless and Gulliver (2018) reported a range of transitional Hf values similar to those in the Nanaimo Basin in Late Cretaceous detrital zircon of the Ochoco Basin and suggested that zircons with transitional Hf values might be derived from the northern Sierra Nevada batholith (Fig. 1). However, transitional Hf values occur in detrital zircon as young as 65 Ma in the Nanaimo Basin, which are too young to have been derived from the Sierra Nevada system. Pecha et al. (2021) documented a distinctive group of Late Cretaceous plutons yielding transitional Hf values (εHf ∼+5 to −2) in the Colorado Mineral Belt that are unique in North America. The Colorado Mineral Belt was the locus of voluminous plutonism, extension, and exhumation in Late Cretaceous time (ca. 75–65 Ma) (Kelley et al., 2001; Chapin, 2012). The Colorado Mineral Belt would provide an ideal source for the transitional Hf zircon in both the Nanaimo and Ochoco basins, although a viable sediment pathway has yet to be identified.

The isotopically evolved character of young zircon of the Nanaimo Basin requires an extra-regional source. The most likely source region capable of producing 75–70 Ma zircon with evolved Hf isotopic values (εHf <−10) is the Atlanta lobe of the Idaho batholith (Fig. 13A; Stroup et al., 2008; Gaschnig et al., 2010). More than 60% of highly evolved young zircon is found in the muscovite-rich Gabriola Formation of the upper Nanaimo Group, consistent with sediment derivation from the Idaho batholith in Maastrichtian–Paleocene time.

Nanaimo Basin Sediment Provenance

The integration of detrital zircon U-Pb geochronology and Lu-Hf analysis, quartzite clast U-Pb geochronology, and muscovite 40Ar/39Ar analysis provides critical constraints on sediment provenance in the Nanaimo Basin. Specifically:

  1. Uranium-lead (U-Pb) geochronology and Hf isotopic constraints on the Mesozoic zircon source: Paleozoic and Mesozoic zircon may be linked to several sources: (a) Devonian–Mississippian zircon were probably derived locally from the Late Devonian Salt Spring Island intrusions and the mid-Paleozoic Sicker Group on the Wrangellia terrane; (b) all pre–80 Ma zircon has a juvenile Hf signature, consistent with derivation from the southern Coast Mountains batholith (Fig. 13A); and (c) 68–80 Ma zircon has three subpopulations, including a juvenile subset attributed to a southern Coast Mountains batholith source, an evolved subset consistent with derivation from the Atlanta lobe of the Idaho batholith, and a transitional subset potentially derived from the Colorado Mineral Belt (Fig. 13A).

  2. Uranium-lead (U-Pb) geochronology of the Proterozoic zircon source: Proterozoic zircon in the Nanaimo Group is consistent with derivation from the Lemhi Group of the Belt-Purcell Supergroup and ca. 1380 Ma plutonic rocks intruding the Lemhi subbasin of central Idaho (Dumitru et al., 2016). Neoproterozoic (ca. 600–500 Ma) zircon is likely derived from Ediacaran to Cryogenian plutons described from central Idaho (Lund et al., 2010). The lack of Grenville-aged zircon suggests that Neoproterozoic strata along the western continental margin did not contribute significant sediment to the Nanaimo Basin. Uplift and erosion of the Lemhi subbasin, ca. 1380 Ma plutonic rocks intruding the subbasin, and Neoproterozoic plutons in Campanian time is supported by identification of detrital zircon derived from these sources in the Late Cretaceous Beaverhead Group of southwestern Montana (Garber et al., 2020).

  3. Quartzite clast U-Pb geochronology: Four out of five quartzite clasts from conglomerates within the Geoffrey Formation are statistical matches for Mesoproterozoic and Cambrian strata in Montana, Idaho, and southern British Columbia (Figs. 10, 11).

  4. Muscovite 40Ar/39Ar analysis: Coarse-grained muscovite from the Gabriola Formation yields 78–76 Ma 40Ar/39Ar detrital ages, consistent with uplift of the source region in late Campanian time (Fig. 12). Late Cretaceous muscovite is unusual in source regions immediately adjacent to the Nanaimo Basin. The Atlanta lobe of the Idaho batholith contains 83–68 Ma two-mica granites rich in both muscovite and zircon and provides an ideal source for voluminous sediment influx into Campanian to Paleocene strata of the Nanaimo Basin.

Analysis of spatial and temporal variations in detrital U-Pb geochronology, 40Ar/39Ar thermochronology, and Lu-Hf isotopic data of the Nanaimo Group demonstrates a bimodal sediment supply that reflects detritus chiefly derived from both the flanking magmatic arc in the southern Coast Mountains batholith and an extra-regional source from the Mesoproterozoic Belt Supergroup and associated plutons as well as the Late Cretaceous Atlanta lobe of the Idaho batholith. These linkages strongly suggest that the Nanaimo Basin was originally located north of the Sierra Nevada–Klamath Mountains system and has subsequently experienced only minimal (<∼750 km) latitudinal translation with respect to the Idaho batholith and Belt Supergroup. However, a definitive resolution of the Baja B.C. translation debate requires examination of other areas along the continental margin that could feasibly yield the same multidimensional provenance signature.

Comparison with the Southwestern Laurentia–Mojave Region

Previous workers have used paleomagnetic, paleobotanic, and provenance data to suggest that the Nanaimo Basin was originally located thousands of kilometers south of its current position (cf. Cowan et al., 1997; Enkin, 2006; Miller et al., 2006; Matthews et al., 2017) at the latitude of southern California or northern Mexico (Enkin, 2006). This reconstruction places the Nanaimo Basin, as well as the combined Intermontane-Insular crustal block, adjacent to the southern end of the Sierra Nevada batholith, the Mojave-Sonora region, and the northern end of the Peninsular Ranges batholith (Fig. 1). If the Nanaimo Basin originated in this southerly region, sediment derived from these prospective source areas should be consistent with provenance data from the Nanaimo Group, including:

  1. Uranium-lead (U-Pb) geochronology and Hf isotopic constraints on the Mesozoic zircon source: Detrital zircon age spectra from Cenomanian to Eocene forearc sequences in southern California are distinct from spectra from the Nanaimo Group (Jacobson et al., 2011). Santonian to Campanian forearc strata in southwestern Laurentia are dominated by 135–85 Ma zircon and are isolated from Mojave-Sonora sources by the mid-Cretaceous arc (Jacobson et al., 2011; Sharman et al., 2015). There is a distinct shift to inboard Mojave sources in Maastrichtian–Paleocene strata, including an influx of latest Cretaceous (85–65 Ma) zircon (Sharman et al., 2015). However, these strata maintain a significant component of mid-Cretaceous (135–85 Ma) zircon. In contrast, Maastrichtian strata of the upper Nanaimo Group are characterized by a unimodal Mesozoic peak of syndepositional zircon (ca. 73 Ma; DZP3, Fig. 8), reflecting rapid exhumation and subsidence.

    Late Cretaceous (ca. 80–72 Ma) plutons are abundant in the Mojave, and latest Cretaceous to Paleocene (ca. 70–50 Ma) granites are found to the east in southern Arizona (southwestern United States) and northern Sonora (northwestern Mexico) (McDowell et al. 2001; Needy et al., 2009; González-León et al., 2011; Chapman, 2017). These Late Cretaceous plutons intrude continental crust and are characterized by evolved Hf isotopic values (Chapman et al., 2018).

    There is a significant difference in the Hf value of Late Cretaceous zircon between southwestern and northwestern Laurentia. Late Cretaceous plutons in the Mojave region yield evolved Hf values (εHf = ∼−3 to −25) (Barth et al., 2016; Howard et al., 2016). Howard et al. (2016) reported the most juvenile Hf value (εHf = ∼+0.5) known from Late Cretaceous plutons in the Mojave region. Conversely, Hf values from Late Cretaceous (ca. 73 Ma average; 81–66 Ma) zircon in the Nanaimo Group form three distinct subpopulations, with the most juvenile group yielding εHf = +11 to +15. There is no known source in the Mojave-Sonora region for Late Cretaceous zircon with the juvenile Hf values documented in the Nanaimo Basin (Fig. 13A).

  2. Uranium-lead (U-Pb) geochronology of the Proterozoic zircon source: Proterozoic zircon from Cenomanian to Eocene forearc sequences in the Mojave region has a distinctly different aggregate detrital zircon signature from that of the Nanaimo Group (Fig. 16). Proterozoic age spectra from both the northern and southern Mojave-Sonora regions display bimodal populations, with a strong peak at ca. 1680 Ma (1750–1630 Ma) and a smaller peak at ca. 1405 Ma (1450–1360 Ma) (Fig. 16). The Nanaimo Group has an essentially inverse bimodal age spectrum, with a major peak at ca. 1380 Ma (1410–1340 Ma) and a smaller peak at ca. 1704 Ma (1430–1675 Ma) (Figs. 9, 16). A cumulative probability plot highlights the dissimilarity among these populations (Fig. 16). Importantly, the strong ca. 1380 Ma peak that characterizes the Proterozoic zircon populations of the Nanaimo Group is essentially nonexistent in the forearc strata from the Mojave region. Analysis of Proterozoic detrital zircon age spectra from Upper Cretaceous forearc strata from the Mojave region indicates that 1400–1320 Ma zircon represent <4% (n= 44/810 grains) of this population. In contrast, 1400–1320 Ma zircon comprises 31% of the aggregate Proterozoic population in the Nanaimo Group.

    Sharman et al. (2015) documented a geomorphic breach in the mid-Cretaceous arc in southwestern Laurentia in Maastrichtian–Eocene time that led to a pronounced transition from local to extra-regional provenance, including an influx of Proterozoic zircon. Prior to this time, Proterozoic zircon was rare in the California forearc. Conversely, the Nanaimo Basin records an extra-regional sediment flux of Proterozoic zircon in early Campanian time, ∼12 m.y. before the provenance shift in southern California.

  3. Quartzite clast U-Pb geochronology: Quartzite cobbles have been identified in Cretaceous strata of the Valle Group in Baja California and from the Nacimiento block in California. Quartzite clasts from the Cenomanian to Santonian Valle Group of the Vizcaino Basin are inferred to be derived locally from Paleozoic miogeoclinal strata or Mesozoic rocks in the southwestern United States (Kimbrough et al., 2006). Quartzite clasts from the Albian to Santonian Nacimiento forearc basin conglomerates are believed to be sourced from Neoproterozoic–Paleozoic shelf strata and Triassic–Jurassic strata of the Colorado Plateau (Johnston et al., 2018). No additional quartzite clast data are available from southwestern Laurentia. Quartzite clasts have also been identified in the Albian to Coniacian Ochoco Basin in central Oregon (Surpless and Gulliver, 2018). Quartzite conglomerate clast detrital zircon populations within the Vizcaino and Ochoco basins are statistically similar and are distinct from those of the quartzite conglomerate clasts analyzed from the Nanaimo Basin (Fig. 11B). There is no obvious source for Nanaimo Group quartzite clasts (n= 4; Fig. 11B) in southwestern Laurentia.

  4. Muscovite 40Ar/39Ar analysis: Muscovite is uncommon in the Sierra Nevada and Peninsular Ranges batholiths, but is common in the Mojave-Sonora region, where it is associated with two-mica granites and underplated schists of the Pelona-Orocopia-Rand schist complex (Jacobson et al., 2011). Late Cretaceous muscovite ages are reported from the schist complex, but this muscovite is fine grained and graphitic and was not exposed until the Miocene (Jacobson et al., 2011).

    Two-mica granite from the Mojave region could potentially yield Late Cretaceous muscovite. However, regional extension in the southwestern Cordillera has resulted in abundant 74–68 Ma muscovite cooling ages in the Mojave region that are inferred to record post-intrusive cooling and exhumation by extensional structures (Wells et al., 2005). Muscovite in the Nanaimo Group predates this regional extension, and, importantly, the sample yielding the oldest muscovite 40Ar/39Ar age (78.33 ± 0.48 Ma) is from the Campanian De Courcy Formation, which predates the geomorphic breach in the southern mid-Cretaceous arc segment and subsequent extra-regional sediment flux into the forearc basin by ∼6–7 m.y.

    In addition, an oft-overlooked consequence of potential linkages between the Nanaimo Basin and the southwestern Laurentia–Mojave region is the implication for Late Cretaceous paleogeography. The Nanaimo Basin accounts for <15% of the entire strike length of the combined Insular-Intermontane block. Restoring the Nanaimo Basin to a position outboard of the modern Mojave region would result in the Insular-Intermontane block overlapping the northern ∼100 km of the Peninsular Ranges batholith and the southern ∼1500 km of the Sierra Nevada batholith. There is no geological evidence for such a paleogeographic configuration in Late Cretaceous time.

Paleogeographic Reconstructions

Comparative analysis of multidisciplinary provenance data from the Nanaimo Basin suggests that source regions in northwestern Laurentia, north of the Sierra Nevada–Klamath Mountains system, are most compatible with available provenance constraints. Proposed tectonic reconstructions that place the Nanaimo Basin, southern Coast Mountains batholith, and Wrangellia south of the Sierra Nevada magmatic system are problematic (Garver and Davidson, 2015; Matthews et al., 2017; Sauer et al., 2017a, 2017b, 2019). However, it is essential to critically examine proposed paleogeographic reconstructions in order to test the validity of previously proposed models.

Proposed Linkages to the Mojave-Sonora Region

Matthews et al. (2017) undertook an extensive detrital zircon analysis of the northern Nanaimo Basin and argued that the most likely source area was the Mojave region of southwestern Laurentia, thus requiring large-magnitude northward translation of the western Canadian Cordillera. Their arguments were based on: (1) the age distribution and timing of influx of Proterozoic zircon populations; (2) an apparent mismatch between Proterozoic zircon of the Belt-Purcell Basin and that of the Nanaimo Basin; and (3) the source of <80 Ma zircon in the Nanaimo Basin.

One of the key arguments linking the Nanaimo Group to the Mojave region is the synchronicity of an influx of Proterozoic zircon in late Campanian–early Maastrichtian time in both the California forearc basin and the Nanaimo Basin. There was indeed a major influx of Proterozoic zircon into the Nanaimo Basin at this time, but this pulse postdates an initial influx of Proterozoic zircon into the southern Nanaimo Basin in the early Campanian. This early Campanian influx of Proterozoic detritus in the Nanaimo Basin was 10–12 m.y. earlier than that in the Mojave region forearc, during a time when the Mojave retroarc was isolated from the forearc by the arc massif (Fig. 7; Jacobson et al., 2011; Coutts et al., 2020). Proterozoic zircon was rare in the Mojave region forearc system prior to Maastrichtian time (Jacobson et al., 2011). In addition, this second pulse of zircon sedimentation in the Nanaimo Basin was accompanied by a major influx of syndepositional zircon (ca. 73 Ma) that is absent in the Mojave region forearc. Importantly, the extra-regional influx of Proterozoic detrital zircon in the Nanaimo Basin was accompanied by a flood of syndepositional Late Cretaceous muscovite, which is also absent in southwestern Laurentia.

The age range of Proterozoic zircon in the California forearc basin and the Nanaimo Basin are similar, but the relative proportion of populations of different ages are quite different (Figs. 16A, 16B). Matthews et al. (2017) argued that it is unlikely that the main Belt-Purcell Basin was the source of Proterozoic zircon in the Nanaimo Group because there should be a higher proportion of Archean zircon and North American magmatic gap (1.60–1.49 Ma) grains in the Nanaimo Group. In addition, they argued that the Nanaimo Group displays a broader range of Mesoproterozoic-age grains (1480–1320 Ma) than the classic “Lemhi doublet” (peaks at ca. 1380 Ma and 1800–1650 Ma) that Dumitru et al. (2016) suggested uniquely fingerprints strata of the Lemhi subbasin of the Belt Supergroup in central Idaho.

The Belt-Purcell Supergroup is not a homogeneous entity but consists of a >20 km stratigraphic succession exposed across >200,000 km2 (Lydon, 2007). The supergroup is structurally and stratigraphically complex and is divided into four groups representing >20 formations in the United States portion alone (Winston and Link, 1993). The different groups tend to yield characteristic detrital zircon age spectra, but these spectra vary laterally and vertically across the basin (Ross and Villeneuve, 2003). The combination of original lateral facies changes and Mesozoic deformation has produced a complex map pattern of the Belt Supergroup that varies dramatically along strike (Winston and Link, 1993). The Lemhi and Missoula groups of the upper Belt-Purcell Supergroup are primarily exposed in the southern portion of the outcrop belt in central Idaho and western Montana, adjacent to the Idaho batholith (Winston and Link, 1993). The Lemhi Group displays an essentially unimodal peak with mode of 1715 Ma (ca. 1760–1670 Ma; n = 1019), a minor peak at 1442 Ma, and a scattering of Paleoproterozoic and Archean grains (Stewart et al., 2010; Link et al., 2016). These rocks are intruded by a series of ca. 1380 Ma anorogenic plutons and sills (Lewis et al., 2010; Gaschnig et al., 2013; Dumitru et al., 2016). Derivative sediments from the Lemhi Group and associated magmatic rocks would be predicted to produce detrital age spectra dominated by ca. 1380 Ma and ca. 1715 Ma peaks. The Proterozoic age spectra of the Nanaimo Group form a bimodal distribution with major peaks at ca. 1380 Ma and ca. 1704 Ma (1720–1700 Ma; n = 320; Figs. 9A, 9B) Matthews et al. (2017) reported major peaks at ca. 1393 Ma and ca. 1702 Ma, with smaller peaks at 1500–1350 Ma and Paleoproterozoic to Archean peaks (Fig. 9), precisely what would be expected from erosion of the Lemhi subbasin.

It is important to note that the doublet of Dumitru et al. (2016) is the result of mixing of Belt-Purcell Supergroup detritus with zircon from slightly younger (ca. 1380 Ma) magmatic rocks, so the strength of the peak and scatter around the peak would be expected to vary. The critical point is the volume of ca. 1380 Ma zircon in the Nanaimo Basin (31% of Proterozoic zircon in this study is 1400–1320 Ma; Matthews et al., 2017, reported 17%). In contrast, Upper Cretaceous strata within the Mojave region forearc basins contain ∼3% ca. 1380 Ma zircon (Jacobson et al., 2011). In addition, the ca. 1380 Ma signal is found in the Nanaimo Basin in strata as old as early Campanian, deposited during a time in which the Mojave region forearc basins were isolated from the Mojave-Sonora region (Sharman et al., 2015).

In addition, the distribution of ca. 1.39 Ga anorogenic granitic plutons in western North America is also important. The majority of these anorogenic plutons in the southwest are older than 1.39 Ga, whereas those in the northwest are younger than 1.39 Ga (Fig. 16C). A comparison of 1350–1450 Ma detrital zircon populations in the Nanaimo Basin and forearc strata from the Mojave region demonstrates that the Nanaimo Basin is dominated by zircon younger than 1.39 Ga, which is consistent with derivation from central Idaho and inconsistent with derivation from southwestern Laurentia (Fig. 16C).

Matthews et al. (2017) argued that <80 Ma zircon grains from the Nanaimo Basin show no evidence of inherited Archean or Proterozoic cores, which are common in the Idaho batholith. However, Boivin (2019) reported that ∼66% of zircons analyzed from Maastrichtian strata in the northern Nanaimo Group contain discrete core-rim domains and that the core dates are predominately Proterozoic (83%) with modes at 1698 Ma and 1388 Ma, which is consistent with core ages reported from the northern Atlanta lobe and the Elk City area of the Idaho batholith (Gaschnig et al., 2013). These cores are interpreted to be derived from the Belt-Purcell Supergroup and associated A-type magmatic rocks (Gaschnig et al., 2013).

Finally, Matthews et al. (2017) argued that widespread Late Cretaceous magmatism in the Mojave region is a good match for <80 Ma zircon populations in the Nanaimo Group. It is true that there was abundant Late Cretaceous magmatism in the Mojave region, but the resulting plutons are characterized by evolved Hf values (εHf = ∼−3 to −25) (Barth et al., 2016; Howard et al., 2016). As previously discussed, Hf values from Late Cretaceous (ca. 73 Ma average; 81–66 Ma) zircon in the Nanaimo Group form three distinct subpopulations (juvenile, transitional, and evolved), with the most juvenile group yielding εHf = +11 to +15. There is no known source for Late Cretaceous zircon with juvenile Hf values in the Mojave-Sonora region (Fig. 13).

Proposed Linkages to the Yakutat Group

Garver and Davidson (2015) argued that the Yakutat Group of the Chugach–Prince William terrane of southeastern Alaska was deposited during Maastrichtian time at the latitude of the Mojave region and outboard of the Nanaimo Basin. They suggested that the eastern part of the Yakutat Group has a unique U-Pb and Hf detrital zircon signature that could only have been derived from a distinctive source in southwestern Laurentia and that it therefore forms a Cordilleran piercing point in Late Cretaceous time, which requires subsequent large-scale latitudinal translation (Garver and Davidson, 2015). The hallmarks of this unique U-Pb and Hf detrital zircon signature are age peaks at ca. 1380 Ma (1388–1365 Ma; n = 9), ca. 1459 Ma (1487–1425 Ma; n = 8; mode recalculated from Garver and Davidson, 2015), and ca. 1722 Ma (1790–1660 Ma; n = 26), all characterized by primarily juvenile Hf values. However, the Lemhi subbasin of central Idaho contains detrital sources with essentially identical age ranges, including ca. 1380 Ma plutonic rocks and the Lemhi Group of the Belt-Purcell Supergroup, which exhibits a bimodal detrital zircon assemblage with a major peak at 1715 Ma (1760–1670 Ma; n = 1019) and a minor peak at 1442 Ma (1475–1410 Ma; n = 68) (Link et al., 2006, 2016; Stewart et al., 2010; Dumitru et al., 2016). Isotopic comparisons are more difficult. There are no known Hf isotopic data for the ca. 1380 Ma grains in central Idaho. There are also no Hf data available for the ca. 1370 Ma anorogenic magmatic zircons in southwestern Laurentia; all analyses of this age are from the Granite-Rhyolite or eastern Yavapai province exposed hundreds of kilometers to the east and not from the Mojave region (Goodge and Vervoort, 2006). The 1460–1440 Ma Hf values are problematic due to limited analyses and overlapping isotopic values, making comparisons difficult. Hafnium values in this age range in the Yakutat Group are juvenile (εHf = +0.4 to +6.8; n = 8), whereas A-type magmatic zircon from the Mojave is both younger (ca. 1420 Ma) and more evolved (εHf = −6 to +5) than the Yakutat population (Wooden et al., 2013). Lemhi Group equivalents in central Idaho overlap the Yakutat Group values but exhibit a broader range of hafnium values (εHf = −8.4 to +4.7; n = 7). Garver and Davidson (2015) pointed out that the 1790–1660 Ma subpopulation of the Yakutat Group is more juvenile (εHf = −3.3 to +9.8) than equivalent zircon in the Lemhi Group of the Belt-Purcell Supergroup (εHf = −7 to +8) and suggested this provides a definitive tie to southwestern Laurentia (Stewart et al., 2010). However, the Hf values of Lemhi Group–equivalent strata vary by stratigraphic position, and the Lawson Creek Formation has Hf values (εHf = +3 to +8) that directly overlap with those of coeval zircon in the Yakutat Group (Link et al., 2016). More importantly, both basement rocks and Mesoproterozoic and Neoproterozoic siliciclastic cover sequences in the Mojave region have Hf values that range between −10 and +10 (Wooden et al., 2013). Hafnium values of Proterozoic zircon from the Yakutat Group therefore overlap with values of potential source areas in both central Idaho and the Mojave region, limiting their utility as unique provenance markers.

There are other aspects of the Yakutat Group detrital zircon age spectra that may constrain potential source terranes. Garver and Davidson (2015) compared the Precambrian grain age distribution of the Yakutat Group with that of the Pelona-Orocopia-Rand schists of southern California and noted that a distinct difference between them is the lack of Grenville-aged (ca. 1200–1000 Ma) zircon in the Yakutat Group. Grenville-aged zircon also forms a very minor component in the Nanaimo Basin and in schists in the North Cascades (Fig. 7; Sauer et al., 2017b, 2018). The Yakutat Group contains a subpopulation of Late Cretaceous (ca. 90–65 Ma) zircon grains that forms three isotopically distinct groups, including juvenile (εHf >+5), transitional (εHf = 0 to +5), and evolved (εHf <−10). This threefold subdivision is similar to that of coeval zircons in the Nanaimo Basin, the Western mélange belt (Fig. 2A, west of North Cascades core), and the Swakane Gneiss (Fig. 13B; Sauer et al., 2018, 2019). We note also that there is no known source for juvenile Late Cretaceous zircon in the Mojave region, whereas this is a distinctive component of the southern Coast Mountains batholith (Fig. 13A; Cecil et al., 2018).

A straightforward interpretation of the origin and displacement of the Yakutat Group is that the group represents a lateral, age-equivalent, basinward facies of the upper Nanaimo Group or other Wrangellia forearc successions (Plafker, 1987; Haggart, 1991b, 1993). The upper Nanaimo Group was deposited on Wrangellia, whereas the eastern Yakutat Group may have been deposited on oceanic crust outboard from Wrangellia that was subsequently displaced to the north by motion along the transcurrent Queen Charlotte fault outboard of the Nanaimo Group in the early Eocene (Cowan, 2003). This interpretation is consistent with tectonic reconstructions of the northward translation of the Yakutat terrane along the Queen Charlotte fault from the latitude of southern Oregon to its current position (Wells et al., 2014; ten Brink et al., 2018).

Proposed Linkages to the Pelona-Orocopia-Rand Schists in Southern California

Sauer et al. (2019) argued that the Swakane Gneiss complex of the North Cascades shares a common sediment source with the Pelona-Orocopia-Rand schists of southern California and suggested that the protolith to the Swakane Gneiss was therefore deposited adjacent to the Mojave in southwestern Laurentia. Proposed provenance ties between the Swakane Gneiss and the Pelona-Orocopia-Rand schists include xenocrystic Proterozoic cores, unradiogenic εHf values for zircons <100 Ma, and a lack of Archean grains. However, the Swakane Gneiss, the Pelona-Orocopia-Rand schists, the Mojave region forearc, and the Nanaimo Basin all contain trace amounts of Archean zircon, so its occurrence is not latitudinally distinctive (Mahoney et al., 1999; Jacobson et al., 2011; Matthews et al., 2017; Sauer et al., 2018, 2019). The wide range of Hf values displayed by 100–81 Ma zircon in the Swakane Gneiss (εHf = +13 to −16) contrasts markedly from those of coeval zircon in the Pelona-Orocopia-Rand schists and requires derivation from a diversity of crustal sources (Fig. 13B). Note that 97% (n= 32/33) of 100–81 Ma zircons in the Pelona-Orocopia-Rand schists are evolved (εHf = −0.4 to −26), whereas 42% (n= 29/69) of coeval zircons in the Swakane Gneiss are juvenile (εHf >0) (Sauer et al., 2017b, 2018, 2019). There is no known source for Late Cretaceous juvenile zircon in the Mojave region. The wide range of Late Cretaceous Hf values is readily available in northwestern Laurentia, including juvenile zircon from the southern Coast Mountains batholith (Homan, 2017) and evolved zircon from the Atlanta lobe of the Idaho batholith (Gaschnig et al., 2011). In addition, the Atlanta lobe and the adjacent Lemhi subbasin yield Late Cretaceous zircon with Proterozoic cores and an abundance of ca. 1380 Ma zircon (Gaschnig et al., 2011; Dumitru et al., 2016). Sauer et al. (2019) also noted that the age peaks and limited range of juvenile εHf values in 160–90 Ma zircon within the Swakane Gneiss match those of the southern Coast Mountains batholith and are not readily available in southwestern Laurentia. It is apparent that there is no detrital signature that uniquely ties the Swakane Gneiss to southern latitudes and that proposed provenance linkages between the Swakane Gneiss and schists in southwestern Laurentia are tenuous at best (Fig. 13B). It is more likely that the Swakane Gneiss represents strata correlative with the upper Nanaimo Group that were underthrust beneath the continental margin in the latest Cretaceous (Sauer et al., 2017a; Fig. 13B).

Preferred Paleogeography: Sediment Provenance North of the Sierra Nevada

Multidisciplinary provenance constraints within the Nanaimo Group document a bimodal detrital signature, with both a locally derived source from the southern Coast Mountains batholith and an extra-regionally derived source from the Atlanta lobe of the Idaho batholith and the Lemhi subbasin of central Idaho (Fig. 17). These provenance constraints are most easily accommodated by tectonic restoration of the Nanaimo Basin, the southern Coast Mountains batholith, and Wrangellia to a position outboard of the Idaho batholith in Late Cretaceous time (Fig. 17). This restoration is consistent with the minimum-fault-offset reconstruction of Wyld et al. (2006). In this mid-Cretaceous (ca. 100 Ma) reconstruction, the southern end of the Insular assemblage restores to the latitude of the Klamath Mountains, and the southwestern Intermontane assemblage straddles the Washington-Oregon border (Wyld et al., 2006). Although the tectonic reconstruction presented in the Wyld et al. (2006) model predates the Nanaimo Basin by ∼10 m.y. and makes assumptions about precise paleogeography, restoration along known fault systems in northwestern Laurentia does indeed place the Nanaimo Basin west of the Idaho batholith in Late Cretaceous time (Fig. 17).

This paleogeographic reconstruction is also consistent with provenance constraints from Cretaceous strata in the Hornbrook-Ochoco Basin, which suggests the Blue Mountains province was located ∼400 km south of its current position, in a position east of the Klamath Mountains and north of the Sierra Nevada in the Late Cretaceous (Surpless, 2015; Surpless and Gulliver, 2018). Gaschnig et al. (2017) supported a southern position for the Hornbrook-Ochoco Basin and argued that U-Pb zircon age spectra and Hf isotope signatures indicate a component of sediment derivation from the Coast Mountains batholith (Fig. 17).

The above paleogeographic reconstruction is consistent with the plate model of Doubrovine and Tarduno (2008), who argued that plate circuit reconstructions, seismic tomography, and mantle flow modeling indicate that the Farallon-Kula Ridge was north of the Sierra Nevada–Klamath Mountains system in Late Cretaceous time. It is important to note that these geologic constraints are also consistent with the moderate-translation hypothesis of Butler et al. (2001), wherein post–mid-Cretaceous latitudinal translation along the Cordilleran margin is limited to ∼1000 km.

This paleogeographic model emphasizes geologic constraints on basin locations, not the paleomagnetic constraints that have been used to infer significant latitudinal offset and which are detailed elsewhere (e.g., Irving et al., 1996; Enkin, 2006). The magnitude of latitudinal displacement is calculated by estimating the difference between paleomagnetic poles estimated from Mesozoic strata in the Insular and Intermontane assemblages and cratonic North America (Fig. 3; Enkin, 2006). Provenance interpretations presented herein establish a robust linkage between the Nanaimo Basin and the Mesoproterozoic Belt Basin of the Intermontane west (Fig. 17). However, if the Belt Basin is itself allochthonous, then large-scale latitudinal displacement of the Nanaimo Basin and its Belt Basin source cannot be disproven (Johnston, 2008; Hildebrand, 2013). Key to the debate is therefore the autochthoneity of the Belt Basin itself. Although it is true that the Belt Basin is displaced within the Cordilleran fold-and-thrust belt, provenance linkages between the Belt Basin and the Laurentian craton have been well established elsewhere. These linkages include: (1) detrital zircon data from the lower two-thirds of the Belt Basin stratigraphy that is consistent with derivation from the western Canadian Shield and the Wyoming craton (Ross and Villeneuve, 2003; Oyer et al., 2014); (2) the presence of 500–490 Ma magmatic zircon grains derived from the Beaverhead magmatic suite in central Idaho in Cambrian strata overlying the Wyoming craton (Link et al., 2017); and (3) the stratigraphic continuity between Ediacaran to Cambrian and Devonian to Mississippian strata overlying the Belt Basin and autochthonous strata overlying the Laurentian craton (McMechan et al., 2021). These linkages demonstrate that the Belt Basin is autochthonous with respect to the Laurentian craton, as is the Nanaimo Basin, which was receiving sediment derived from the Belt Basin in Late Cretaceous time.

With respect to the broader issue of the validity of the Baja B.C. hypothesis, it is not within our purview to explain the reason for the discrepancy between the geologic and paleomagnetic data sets. However, the large-scale (>2000 km) latitudinal displacements inferred from paleomagnetic data from the Canadian Cordillera are incompatible with the multidisciplinary geological constraints on the evolution of the Nanaimo Basin documented herein.

Basin Evolution Synthesis

The Nanaimo Group provides an excellent record of orogenic exhumation and basin evolution and constrains the tectonic evolution of the combined Insular and Intermontane crustal blocks (Figs. 17). Subsidence in the Nanaimo Basin began in the early Turonian (ca. 92–90 Ma), coincident with a major period of west- and southwest-vergent contraction and associated magmatism in the Coast Mountains batholith (Journeay and Friedman, 1993; Umhoefer and Miller, 1996). Initial basin sedimentation onlapped a complex paleotopography, resulting in accumulation of locally derived sediment in nonmarine, marginal-marine, and shallow-marine environments. Regional contraction uplifted Middle Jurassic–Early Cretaceous plutons and pre- and syn-kinematic plutons in the southern Coast Mountains batholith that are consistent with the broad Jurassic–Cretaceous trimodal detrital zircon population evident in DZP1 of the lower Nanaimo Group succession (Figs. 7, 8). Locally derived Devonian–Mississippian sediment may have resulted from extensional(?) uplift of the outer arc high in Turonian–Santonian time (ca. 93–85 Ma) (England and Bustin, 1998).

There was a major tectonic reorganization in the southern Coast Mountains batholith and North Cascades in Late Santonian–early Campanian time, characterized by the transition from significant crustal thickening, doubly vergent contraction, plutonism, and peak metamorphic conditions between ca. 96–84 Ma to the initiation of dextral transpression and exhumation at ca. 84–82 Ma (Rusmore and Woodsworth, 1991; Journeay and Friedman, 1993; Umhoefer and Miller, 1996; Brown et al., 2000; McClelland et al., 2000; Brown, 2012; Mitrovic, 2013; Monger and Gibson, 2019; Coutts et al., 2020).

The transition to dextral transpression and exhumation in late Santonian–early Campanian time was coincident with a significant increase in subsidence in the Nanaimo Basin and progradation of deep-marine submarine-fan complexes into the basin (Mustard, 1994; Mahoney et al., 2003, 2016; Englert et al., 2020; Coutts et al., 2020). The submarine channel system appears to have developed several million years earlier in the south, as documented by an influx of channelized coarse clastic strata in the Protection Formation (Mustard, 1994; Mahoney et al., 2016; Coutts et al., 2020). The sedimentation rate during Campanian time, calculated using stratigraphic thicknesses within the northern Nanaimo Basin, was ∼0.10 mm/yr, although this must be considered a minimum due to stratigraphic architecture that suggests sediment bypass was an important process during this interval (Figs. 14, 15; Englert et al., 2020).

The rapid increase in subsidence and onset of deep-marine sedimentation in early Campanian time was coincident with a major shift in sediment provenance. A pulse of exhumation in the southern Coast Mountains batholith exposed both syn- and post-kinematic plutons, yielding the bimodal ca. 88 Ma and ca. 75 Ma detrital zircon signature evident in DZP2 (Figs. 7, 8). This exhumation was accompanied by a major breach in the arc massif and delivery of extra-regional detritus into the southern Nanaimo Basin. Significant (locally >50%) quantities of Proterozoic detrital zircon (ca. 1800–1300 Ma), derived from the Belt-Purcell Supergroup and associated magmatic rocks to the east, mixed with Jurassic–Cretaceous grains derived from the Coast Mountains batholith (Figs. 7, 8). This flux of older zircon corresponds to late Santonian–early Campanian uplift of Mesoproterozoic successions along basement-involved thrust sheets of the Sevier orogenic system to the east (DeCelles, 2004; Fuentes et al., 2011; Laskowski et al., 2013; Garber et al., 2020).

There was a significant shift in basin dynamics in late Campanian–early Maastrichtian time, as evidenced by a major increase in sedimentation rate and a distinct shift in sediment provenance (Figs. 14, 15). The late Campanian–early Maastrichtian (ca. 74–70 Ma) time interval corresponds to maximum crustal thickening in the central Coast Mountains batholith (Mahoney et al., 2009), the culmination of a high-flux magmatic event in the southern batholith (Cecil et al., 2018), initiation of exhumation of the Central Gneiss Complex in the central batholith (Rusmore et al., 2005), and maximum burial and metamorphism of underthrust forearc sediments in the North Cascades (Matzel et al., 2004; Sauer et al., 2017b, 2018). Importantly, inboard from the main Cretaceous arc system, major displacement along the Lewis-Eldorado megathrust marks the transition between Sevier thin-skinned–style deformation and Laramide thick-skinned–style deformation, which is thought to reflect increased coupling between the North American and Farallon plates during flat-slab subduction (Yonkee and Weil, 2015). This sharp increase in regional tectonism occurred during a reorientation of plate geometry that resulted in more oblique convergence and higher tangential velocity of the Kula plate (Doubrovine and Tarduno, 2008; Coutts et al., 2020). The synchroneity of crustal thickening and magmatic episodicity within the arc, hinterland tectonism, and forearc basin subsidence suggest a direct linkage between subduction dynamics and upper-crustal response.

This dramatic increase in regional tectonism resulted in a near doubling of the sedimentation rate in the Nanaimo Basin (0.10–0.19 mm/yr; Fig. 14), a major pulse of coarse-clastic sedimentation represented by pebble to boulder conglomerate of the Geoffrey Formation, a shift from lateral channel migration to vertical aggradation in coarse-clastic channelized facies, and the development of multiple fan-channel networks and associated submarine conduits sourced in topographically elevated areas in both the arc and the hinterland (Mahoney et al., 2016; Bain and Hubbard, 2016; Englert et al., 2018, 2020). The increased sediment flux was characterized by a distinct provenance shift among both Mesozoic and Proterozoic populations (Fig. 16). The Mesozoic population is characterized by an abrupt shift to an essentially unimodal peak of 75–68 Ma detrital zircon in the Northumberland, Geoffrey, and lower Gabriola (central transect) formations, a peak that becomes only slightly more subdued in the Spray and Gabriola formations (Figs. 7, 8). This flood of syndepositional grains is important because the very short lag time (∼2–3 m.y.) between crystallization and deposition of the 75–68 Ma zircon requires rapid exhumation of the source plutons. In addition, Hf isotopic data from these grains requires derivation from at least two distinct source areas, including both the southern Coast Mountains batholith and the Atlanta lobe of the Idaho batholith, indicating rapid exhumation in both the arc and the hinterland in late Campanian–early Maastrichtian time. The Proterozoic population increased to ∼26% of the total detrital zircon population (Figs. 8, 9), and, importantly, the ratio of ca. 1380 Ma grains to 1800–1650 Ma grains increased dramatically, from ∼1:1 in strata between late Santonian–late Campanian to ∼2:1 in strata of late Campanian–Maastrichtian age. This increase in ca. 1380 Ma grains, together with the influx of extra-regional quartzite clasts in conglomeratic facies, is consistent with major uplift and displacement along the Lewis-Eldorado megathrust system in western Montana, which exhumed the Belt-Purcell Supergroup, including the Lemhi subbasin, and associated magmatic rocks (Fuentes et al., 2011; Dumitru et al., 2016).

The late-stage breaching of the arc massif and subsequent delivery of extra-regional detritus into the depocenter evident in the Nanaimo Basin is a characteristic attribute of Cordilleran forearc basin evolution. A similar geomorphic breach is well documented in southwestern Laurentia, where there is a pronounced transition from local to extra-regional provenance, including an influx of Proterozoic zircon into the forearc in Maastrichtian–Eocene time (Sharman et al., 2015). Similarly, Finzel et al. (2016) documented a transition from proximal arc sources to sources in the retroarc in southern Alaska during the early Paleogene. These local to extra-regional provenance shifts in each area reflect the progressive increase in surface elevation of backarc regions between Late Cretaceous and middle Eocene time, coincident with Sevier-Laramide deformation. The building and demise of this orogenic highland was thus a continental-scale event extending many thousands of kilometers along strike. Multiple geodynamic models have been proposed as the cause of this surface uplift: upper-crustal response to transitory reorganizations in subduction dynamics (Nanaimo Basin), interaction of a subducted oceanic plateau with continental crust (Mojave region), and subduction of a spreading ridge (southern Alaska). Provenance studies of forearc basin successions, including detrital zircon geochronology and low-temperature thermochronology, provide key data in documenting changing geodynamic conditions and upper-crustal response within arc systems (Finzel et al., 2016; Fasulo et al., 2020; Isava et al., 2021).

The Turonian to Paleocene Nanaimo Group provides a high-fidelity record of orogenic exhumation and basin evolution in the southern Canadian Cordillera (Fig. 1). A comprehensive synthesis of stratigraphic architecture, biostratigraphy, and detrital zircon geochronology facilitates a detailed reconstruction of the Nanaimo Group depocenter and constrains the linkage between episodic basin subsidence, structural deformation, and orogenic exhumation patterns in the southern Canadian Cordillera.

This investigation integrates detrital zircon U-Pb geochronology, conglomerate clast U-Pb geochronology, detrital muscovite 40Ar/39Ar thermochronology, and zircon Lu-Hf isotopic analysis of detrital zircon to constrain the source region(s) supplying sediment to the Nanaimo Basin in Late Cretaceous time. These multidisciplinary provenance constraints require influx from multiple sources, including a locally derived source from the southern Coast Mountains batholith and an extra-regionally derived source from the Atlanta lobe and Lemhi subbasin of central Idaho. These provenance constraints are most easily accommodated by tectonic restoration of the Nanaimo Basin, the southern Coast Mountains batholith, and Wrangellia to a position outboard of the Idaho batholith in Late Cretaceous time.

The results of this investigation are clear: The Nanaimo Basin is a Late Cretaceous forearc depocenter that accumulated north of the Sierra Nevada–Klamath Mountains system, west/northwest of the Idaho-Boulder batholith and a few hundred kilometers southeast of its current geographic position. At this point, the outstanding issue with the Baja B.C. hypothesis is how to explain the consistent paleomagnetic data that are interpreted to indicate large-scale latitudinal displacement.

This manuscript is the result of a multi-year, international interdisciplinary investigation involving several universities and the Geological Survey of Canada. Peter Mustard (Simon Fraser University, retired) provided the original inspiration for the investigation and was responsible for the majority of the stratigraphic analyses. The majority of this investigation was conducted between 1998 and 2010, assisted by a number of undergraduate students from the University of Wisconsin–Eau Claire, including C. Rowe, J. Morrison, M. Haskin, A. Wiest, and others. Kirsten Hodge (University of British Columbia) and Betty Franklin and Dan Bowen (Courtenay, British Columbia) are thanked for help with sample collection. Carl Jacobson is thanked for his contribution and insight into the tectonic evolution of southwestern Laurentia. The study was funded by U.S. National Science Foundation grants EAR-9628515 and EAR-0409869 and by support from the Office of Research and Sponsored Programs at the University of Wisconsin–Eau Claire. J.W.H Monger, L. Beranek, and J. Nelson are thanked for discussions and insights that clarified and enhanced our presentation herein. Jeff Trop, Ken Ridgway, and Steve Johnston provided thorough reviews and significantly improved the manuscript. This publication is long overdue, which is the sole fault of the senior author, who is also responsible for any errors or omissions in the manuscript.

1Supplemental Material. Item A: Summary of Nanaimo Group Detrital Zircon Samples (A1–A2). Item B: Biostratigraphic Age Control on Samples (B1–B10). Item C: Nanaimo Group Detrital Zircon Data (C1–C77). Item D: Nanaimo Group Conglomerate Clast Detrital Zircon Data (D1–D7). Item E: Nanaimo Group Muscovite Ar-Ar Geochronology (E1–E7). Item F: Nanaimo Group Lu-Hf Isotopic Data (F1–F7). Item G: Modern Salmon River Detrital Zircon Data (G1–G3). Please visit https://doi.org/10.1130/GEOS.S.16455609 to access the supplemental material, and contact editing@geosociety.org with any questions.
Science Editor: Andrea Hampel
Associate Editor: Christopher J. Spencer
Gold Open Access: This paper is published under the terms of the CC-BY-NC license.