Jurassic–Cretaceous sedimentary rocks of the Alberta foreland basin are a key record of the evolution of the Canadian Cordillera. We test a recent model for cyclical development of Cordilleran orogenic systems using detrital zircon analysis of the major sandstone units deposited between 145 and 80 Ma exposed in the Rocky Mountain Foothills near Grande Cache, Alberta. The basin history is well constrained by decades of study, and the stratigraphy has been previously subdivided into tectonostratigraphic wedges. U-Pb data from 14 detrital zircon samples are included in this study. All the major magmatic provinces of North America are represented in each sample, with the relative proportions varying between samples. The samples are assigned to five groups with the aid of multidimensional scaling. Groups 1–3 are interpreted to record recycling from specific passive-margin units of western North America with varying input from the Cordilleran magmatic arc. Group 4 is interpreted to record recycling from sedimentary strata in the United States and dispersal by basin-axial fluvial systems. Group 5 is dominated by Mesozoic zircon grains interpreted to have originated in the Cordilleran magmatic arc. Detrital zircon age spectra do not form groups based on the tectonostratigraphic wedges from which they were sampled; rather, within each tectonostratigraphic wedge, they exhibit evolution from diverse age spectra to a less-diverse distribution of detrital zircon ages.

We constructed a proxy for magmatic flux of the Cordilleran magmatic arc using detrital zircon ages younger than 200 Ma; it shows three modes at ca. 165, 115, and 74 Ma. These ages are considered high-flux episodes of magmatism that are linked to cyclical uplift and plateau formation in the orogen. This cyclical process is interpreted to: (1) control sedimentation rates in the foreland; (2) account for evolving provenance by altering catchments; and (3) be a plausible mechanism for the deposition of the tectonostratigraphic wedges in the Alberta foreland basin.

Recent models of cyclical orogenic development attempt to describe a unifying framework of interrelated crustal processes, including underthrusting, eclogite root foundering, crustal shortening, episodic magmatism, and plateau development and collapse (DeCelles et al., 2009; Vanderhaeghe, 2012). Inherently, these models predict episodic sedimentation in foreland basins. As such, foreland basin strata are a key archive in which to sample Cordilleran magmatic arcs and test orogenic cyclicity hypotheses.

The western North American foreland basin extends from southern Mexico to the Canadian Arctic, ∼6000 km along strike, with a maximum width exceeding 1000 km (DeCelles, 2004). The Alberta foreland basin is the portion of this basin occupying the Canadian province of Alberta. The Canadian Cordillera and the Alberta foreland basin are exceptionally well studied due to expansive outcrops and hundreds of thousands of well penetrations. This linked orogen-basin system is the focus of classic works on accretionary margin tectonics, fold-and-thrust belts, basin analysis, stratigraphy, and sedimentology (Bally et al., 1966; Monger et al., 1972; Monger and Price, 1979; Coney et al., 1980; Beaumont, 1981; Monger et al., 1982; Stott, 1984; Cant and Stockmal, 1989; Leckie and Smith, 1992).

The Alberta foreland basin formed and filled in response to tectonic loading of the western margin of North America by allochthonous and parautochthonous terranes starting in the Middle Jurassic (Monger et al., 1972, 1982; Monger and Price, 1979). Docking of terranes to the North American margin progressed until the Eocene, resulting in a complex orogenic collage (Monger et al., 1972; Coney et al., 1980).

Numerous authors have recognized the cyclical nature of Alberta foreland basin strata and subdivided the fill into lithostratigraphic cycles or tectonostratigraphic wedges (Stott, 1984; Cant and Stockmal, 1989; Leckie and Smith, 1992; Ross et al., 2005; Pana and van der Pluijm, 2015). The paleogeographic evolution of the basin is well understood, providing constraints on accommodation development and sediment-routing variation (Jackson, 1984; Leckie and Smith, 1992).

Despite this well-established framework for the Alberta foreland basin, the origin of sediments throughout the stratigraphic section remains poorly documented. In this study, we analyzed every significant sandstone unit deposited between 145 and 80 Ma from the Grande Cache area of the central Alberta Foothills using U-Pb geochronology of detrital zircon grains. We tested the hypothesis of orogenic cyclicity in western Canada through analysis of detrital zircon spectra in the context of a high-resolution stratigraphic framework from this uniquely well-constrained foreland basin.

The stratigraphic framework for siliciclastic Mesozoic units in the Alberta foreland basin has been extensively analyzed, with several studies emphasizing the linkage of sedimentary packages to tectonic processes in the adjacent Canadian Cordillera (Fig. 1; Table 1; Cant and Stockmal, 1989; Ross et al., 2005; Raines et al., 2013; Pana and van der Pluijm, 2015). The stratigraphy of the Alberta foreland basin consists of unconformity- or flooding surface–bounded sequences, which are variably subdivided into tectonostratigraphic wedges or lithostratigraphic cycles (Fig. 1; Stott, 1984; Cant and Stockmal, 1989; Leckie and Smith, 1992; Ross et al., 2005; Pana and van der Pluijm, 2015). These cycles provide the framework for the timing of sedimentation and hiatus events, which can be compared to tectonic events in the Canadian Cordillera.

Jurassic to earliest Cretaceous deposits in the basin are widely assigned to the first tectonostratigraphic wedge or depositional cycle in the basin. Protracted subsidence and sedimentation are linked to loading of the lithosphere by accretion of the Intermontane superterrane to the western margin of North America (Cant and Stockmal, 1989). In the study area, these deposits are represented by the Monteith Formation of the Minnes Group, which represents a progradational package of deltaic to fluvial sediments (Miles et al., 2012; Kukulski et al., 2013a).

The basal bounding surface of the second tectonostratigraphic wedge is the basinwide sub-Cretaceous unconformity, which represents a 10–20 m.y. hiatus attributed to isostatic rebound during an extended period of tectonic quiescence (Heller et al., 1988; Cant and Stockmal, 1989). The Bullhead Group, basal Fort St. John Group, and equivalents are assigned to the second cycle of sedimentation in the basin (Table 1; Cant and Stockmal, 1989; Leckie and Smith, 1992; Ross et al., 2005). Accretion of the Bridge River terrane to the margin of North America at this time has been linked to deposition of these sediments (Price et al., 1981; Rusmore et al., 1988; Cant and Stockmal, 1989).

There is disagreement as to whether the Cadotte Member of the Peace River Formation (Fort St. John Group) is part of the second tectonostratigraphic wedge or should be considered as part of an intervening period of tectonic quiescence (Table 1; Cant and Stockmal, 1989; Leckie and Smith, 1992). Gouge of comparable age was absent from major thrust faults in the Rocky Mountains, consistent with the tectonic quiescence hypothesis (Pana and van der Pluijm, 2015).

The Dunvegan Formation represents a major deltaic system that is assigned to the third tectonostratigraphic wedge (Fig. 1; Cant and Stockmal, 1989; Bhattacharya, 1994; Plint, 2000). Deposition of the Dunvegan Formation has been linked to accretion of the Cascadia terrane and a period of oblique compression (Price et al., 1981; Cant and Stockmal, 1989; Pana and van der Pluijm).

The Smoky Group is not included in the third tectonostratigraphic wedge but is rather linked to a time of quiescence in the orogen (Cant and Stockmal, 1989). There are no fault gouge ages that overlap this time interval (Pana and van der Pluijm, 2015).

A fourth major clastic wedge is widely reported for the basin (Fig. 1). Campanian to Paleocene sediments of the Saunders Group and equivalents, including the Belly River Group, Edmonton Group, and Paskapoo Formation, record the last major phase of foreland sedimentation (Cant and Stockmal, 1989; Leckie and Smith, 1992; Ross et al., 2005). This stratigraphy is dominantly nonmarine and has been linked to accretion of the Insular superterrane, Pacific Rim terrane, and Olympic terrane (Price et al., 1981; Cant and Stockmal, 1989; Leckie and Smith, 1992). This tectonostratigraphic wedge was not sampled because the strata are not exposed in the Grande Cache area.

Grande Cache, located in the Rocky Mountain Foothills of west-central Alberta, is an ideal place in which to evaluate the provenance evolution in the Alberta foreland basin for three reasons: (1) Nearly every significant sandstone of the foreland basin sequence from 145 to 80 Ma crops out within 15 km of the townsite. (2) Research focused on the structural geology of the area constrains the distribution and geological history of the units investigated (e.g., Mountjoy, 1980; McMechan, 1999). (3) Various stratigraphic, sedimentologic, and paleontologic studies have characterized units in the area (Fig. 2; e.g., Plint, 2000; Collom, 2001; Kukulski et al., 2013b; McCrea et al., 2014).

Strata were sampled within 15 km of Grande Cache, primarily within the Muskeg thrust (Fig. 2; Irish, 1951; Richardson et al., 1991; McMechan, 1996; McMechan and Wozniak, 1996). The stratigraphic context for samples analyzed is presented in Figure 3.

Four sandstones of Jurassic to earliest Cretaceous age from the Monteith Formation of the Minnes Group are included in this study. The informal Monteith C and Monteith B members were analyzed as part of this study (Miles et al., 2012). These data are augmented with previously published data from the Monteith C and A members (Raines et al., 2013). Two samples were analyzed from the Barremian–Aptian Bullhead Group, including the Cadomin and Gething Formations. Five samples were collected from Albian–Cenomanian sandstones of the Fort St. John Group, including the Bluesky Formation, the Falher and Notikewin Members of the Spirit River Formation, and the Cadotte Member of the Peace River Formation. One sample was collected from the Cenomanian Dunvegan Formation and two from the Cenomanian–Campanian Smoky Group (Cardium Formation and the Chinook Member of the Puskwaskau Formation; Fig. 3).

Sandstone samples, 3 to 4 kg each, were pulverized, and density separation was performed on an MD Gemini Goldharvester™ shaking table (water table). Zircon grains were further concentrated using heavy liquids. Magnetic separation was followed by dump-mounting of grains in 2.5-cm-diameter epoxy puck mounts, which were ground to expose the zircon grains and polished. U-Pb age analysis was performed at the Centre for Pure and Applied Thermochronology and Tectonics (CPATT) at the University of Calgary using laser ablation–inductively coupled plasma–mass spectrometry (LA-ICP-MS). The zircon grains were analyzed using a 22 µm spot diameter. Four zircon references were ablated along with the unknowns to correct for laser-induced elemental fractionation, instrumental fractionation, drift, and to test the accuracy of the laser-ablation procedure (Temora-2—Black et al., 2004; 91500—Wiedenbeck et al., 1995; FC-1—Paces and Miller, 1993; 1242—Mortensen and Card, 1993). Data reduction was performed in Iolite™ V2.5 following the methodology of Paton et al. (2010), using components of the Vizualage data reduction scheme (Petrus and Kamber, 2012). Final assessment of uncertainty and data visualization were done in Excel™ using Isoplot (Ludwig, 2012). All ages are reported as concordia ages (Ludwig, 1998). Errors are quoted at the 2σ level. Full isotopic data and details of the analytical method can be found in the GSA Data Repository.1

Maximum depositional ages were calculated as the weighted average of the youngest subset of zircon grains that overlapped in age at 2σ error (Dickinson and Gehrels, 2009). At minimum, three overlapping zircon ages were used in each case. The youngest single grain in each sample is also reported, with the corresponding error at the 2σ level. Although this method has been widely applied, it perhaps incorrectly assumes that the youngest overlapping grains are comagmatic (cf. Spencer et al., 2015a). Although future analyses should perhaps consider this approach, the method of Dickinson and Gehrels (2009) is deemed adequate for the level of interpretation utilized in this study.

Multidimensional scaling (MDS) is used to aid in the visualization and interpretation of statistical differences in detrital zircon spectra. MDS creates a map of the data points based on the D values that are produced as part of the Kolmogorov-Smirnov test (Vermeesch, 2013; Spencer and Kirkland, 2016). The D value represents the maximum vertical separation between two cumulative density functions. These “distances” are arranged in a matrix, which is then represented in two-dimensional space such that all the “distances” in the matrix are honored. Uncertainties are not considered in this analysis. Limitations of using the D value include the following: (1) It is known to have higher sensitivity to differences in the cumulative density functions near the median; (2) it is more sensitive to proportions of ages rather than the mean age of the mode; and (3) it is not as robust when comparing small data sets (n < 300–600; Saylor and Sundell, 2016). Despite these limitations, MDS provides reasonable groupings in this data set. Unimodal, normally distributed age populations were created in Excel™ and added to the MDS map to show key inputs to the various groups. A Matlab™ code and a graphical-user interface made available by Vermeesch (2013) were used to perform the MDS analysis and create the plots. The matrix of D values from the Kolmogorov-Smirnov test is included in the data repository item (see footnote 1).

Prominent detrital zircon age modes, maximum depositional ages, depositional environments, paleogeographic context, tectonic significance, and concurrent orogenic events are reviewed for each of the stratigraphic horizons sampled in Table 1. The age modes referenced in Table 1 are reported in order of prominence. Probability density plots for each of the samples are presented in Figure 4.

Jurassic–earliest Cretaceous strata that comprise the first tectonostratigraphic wedge have age modes spanning the Mesozoic to Archean (Fig. 4). Similarities between the samples include modes in the Paleozoic (ca. 430 Ma), early Mesoproterozoic (ca. 1000 Ma), the Paleoproterozoic (ca. 1800 Ma), and the Archean (ca. 2700 Ma). The major difference between these samples is the proportion of ages from the various age groups. The Monteith C member analysis of Raines et al. (2013) is unique in that the Proterozoic is characterized by more numerous and prominent modes than the other samples. This is evident in important modes centered at ca. 1050 Ma, ca. 1325 Ma, ca. 1510 Ma, and ca. 1790 Ma. The Monteith C and Monteith B member samples are unique in that Paleozoic and Neoproterozoic modes are more important constituents, including ages at ca. 550–600 Ma. The Monteith A member sample analysis of Raines et al. (2013) is unique in that it is dominated by a strong age mode in the Paleoproterozoic at ca. 1840 Ma (Fig. 4).

Early Cretaceous strata comprising the second tectonostratigraphic wedge have detrital zircon ages spanning the Mesozoic to the Archean (Fig. 4). All six samples share a similar age mode in the Jurassic at ca. 160 Ma. Excluding the Notikewin Member sample, age modes in the Paleozoic (ca. 400 Ma), the early Mesoproterozoic (ca. 1000 Ma), and the Paleoproterozoic (ca. 1800 Ma) are also common between the samples. The Cadomin Formation sample is dominated by a Cretaceous mode at 117 Ma and a Paleoproterozoic mode at ca. 1840 Ma. The Gething Formation has similar modes to the Cadomin Formation at 117 Ma and ca. 1790 Ma; however, the Gething Formation has a more prominent early Mesoproterozoic mode at 1108 Ma. The Bluesky Formation is most similar to the Cadomin Formation in that it is dominated by modes in the Cretaceous (114 Ma) and the Paleoproterozoic (ca. 1760 Ma and ca. 1800 Ma). The overlying Falher D Member has a slightly broader spectrum of ages in the Proterozoic and a higher proportion of Mesoproterozoic zircon grains than the Bluesky Formation. The Falher A Member is unique due to the high proportion of Paleozoic to Neoproterozoic grains (300–600 Ma). Likewise, the Mesoproterozoic mode at 1040 Ma is more prominent in the Falher A Member than the Paleoproterozoic mode at ca. 1850 Ma. The Notikewin Member is dominated by Mesozoic zircon grains with modes at 165 Ma, ca. 110 Ma, and ca. 210 Ma. Proterozoic detrital zircon grains with similar ages to the other samples were analyzed from the Notikewin Member sample, but the modes are muted because of the number of Mesozoic zircon grains (Fig. 4).

Early–Late Cretaceous strata include age modes spanning the Mesozoic to the Archean (Fig. 4). The four samples also share age modes in the Mesozoic (110–120 Ma), the Paleozoic (400–500 Ma), the Mesoproterozoic (1000 Ma), and the Paleoproterozoic (1600–1800 Ma). The Cadotte Member and Dunvegan Formation contain lesser proportions of Mesozoic zircon grains than the Cardium Formation and Chinook Member samples. The Cadotte Formation exhibits a prominent mode in the Paleozoic at ca. 440 Ma, and the Proterozoic zircon grains are evenly distributed between the Mesoproterozoic and Paleoproterozoic, with age modes at ca. 1670, ca. 1200, ca. 1730, ca. 1050, and ca. 1490 Ma, in order of decreasing prominence. The Dunvegan Formation is distinct from the Cadotte Member in that the mode at ca. 1770 Ma dominates the Proterozoic detrital zircon spectra. The Proterozoic detrital zircon grains in the Chinook Member and Cardium Formation are evenly distributed between the Mesoproterozoic and the Paleoproterozoic, though the mode at 1800 Ma in the Cardium Formation is somewhat more prominent (Fig. 4).

Incorporation of progressively younger zircon grains occurs up stratigraphic section (Fig. 4). However, the youngest age mode is not always dominant. In the case of the early Albian Notikewin Member, Middle Jurassic detrital zircon grains compose most of the sample; the Cenomanian Cardium Formation has a most prominent Aptian age mode (120 Ma). In the Chinook Member, however, near-depositional-age detrital zircon grains dominate (84 Ma).

The ultimate source areas of detrital zircon grains in these strata reflect the Archean–Mesozoic magmatic history of North America. Archean zircon grains originated during the formation of the cratonic nucleus of North America (Hoffman, 1988). Paleoproterozoic detrital zircon grains are associated with several orogenic sources that sutured the Archean cratons of North America (Hoffman, 1988; Whitmeyer and Karlstrom, 2007). Detrital zircon age modes between 1800 and 1900 Ma are referred to as Trans-Hudson grains, though there are other Paleoproterozoic orogens that have magmatic assemblages of a similar age (Hoffman, 1988; Whitmeyer and Karlstrom, 2007). The Yavapai-Mazatzal Provinces contributed Paleoproterozoic zircon grains that are younger than the Trans-Hudson orogen, at 1610–1790 Ma (Hoffman, 1988). Mesoproterozoic detrital zircon grains with ages from 1340 to 1480 Ma originated from the midcontinent belt of anorogenic plutons known as the Granite-Rhyolite Province (Windley, 1989). The Grenville orogen (ca. 1000–1300 Ma) contributed Mesoproterozoic detrital zircon grains that are common in this study; the Grenville magmatic pulse at 1040 Ma is particularly well represented (Easton, 1986; Dickinson, 2008). Younger Proterozoic zircon ages are associated with allochthonous terranes of the Appalachian orogeny, specifically the Brasiliano and Pan-African orogens of Gondwana (550–850 Ma; Park et al., 2010). Other Appalachian and pre-Appalachian rift magmatic assemblages span the Paleozoic to Neoproterozoic (285–760 Ma; Dickinson and Gehrels, 2009; Park et al., 2010).

Mesozoic detrital zircon age modes can be attributed to the magmatic arcs and accreted terranes that constitute the Cordilleran orogen to the west. Triassic plutons do not exist on pre-Cordilleran North America; however, igneous and sedimentary units of this age are common to accreted terranes (e.g., Armstrong, 1988; LaMaskin, 2012). Major phases of arc magmatism in the Jurassic and Early Cretaceous occurred at 165–170 Ma and ca. 100 Ma in the Omineca belt of interior British Columbia (Archibald et al., 1983). Late Cretaceous magmatism is represented in the Idaho batholith, the Shuswap complex of southeastern British Columbia, and isolated areas of central British Columbia south of the Bowser Basin (Hyndman, 1983; Armstrong, 1988; Carr, 1995).

In most instances, the biggest difference between individual zircon age spectra is not the presence of a particular mode but the proportions of age groups (Fig. 4). Many of the age spectra have relatively evenly distributed proportions among the Mesozoic, Paleozoic, and the three eras of the Proterozoic. There are important exceptions, however. For example, the Monteith A member is composed of 51% Trans-Hudson–age detrital zircon grains, and the Notikewin Member is composed of 75% Cordilleran zircon grains. The Cardium Formation and Chinook Member are less dominated by Cordilleran influence, composed of 16% and 18% Cordilleran zircon grains, respectively.

The detrital zircon age spectra do not form groups based on their assignment into lithostratigraphic groups nor the tectonostratigraphic wedges that have been interpreted by numerous authors (Figs. 1 and 3; Stott, 1984; Cant and Stockmal, 1989; Leckie and Smith, 1992; Ross et al., 2005). Jurassic to earliest Cretaceous rocks of the first tectonostratigraphic wedge are characterized by detrital zircon grains that evolve from a broad age spectra to one dominated by Trans-Hudson ages (Fig. 4). Following the 10–20 m.y. hiatus of the sub-Cretaceous unconformity, the second tectonostratigraphic wedge evolves from broad age distributions in the Bullhead Group and basal Fort St. John Group to a Mesozoic-dominated signature in the Notikewin Member (Fig. 4). The units from the uppermost stratigraphy studied exhibit a similar pattern, evolving from broad age spectra in the Cadotte Member and Dunvegan Formation to one dominated by Mesozoic modes (Fig. 4).

In the context of the assembly of North America, the Cordilleran orogen and foreland basin system are young. The ultimate zircon sources can be linked to the Alberta foreland basin directly via continent-scale Jurassic–Cretaceous–age sediment routing systems (e.g., Raines et al., 2013; Benyon et al., 2014, 2016; Blum and Pecha, 2014) or indirectly through one or more phases of sediment recycling (e.g., Dickinson and Gehrels, 2009; Gehrels and Pecha, 2014; Hadlari et al., 2015). Recycling of sedimentary rocks imparts the full detrital zircon age spectrum from the eroded strata to the subsequent deposit; therefore, comparison of detrital zircon spectra from older and younger sedimentary rocks is a valid means of detecting sedimentary recycling (Hadlari et al., 2015). As a result, consideration of continentwide, pre-Jurassic sediment routing and recycling is critical for interpreting the provenance of detrital zircon grains in the study area.

The western margin of North America is host to thick packages of sedimentary rocks that accumulated in rift, passive-margin, and foreland basin settings from the Mesoproterozoic to the Triassic (Price et al., 1972; Root, 2001; Dickinson, 2004). These strata were uplifted during the Cordilleran orogeny and contributed a diverse array of detrital zircon grains to the Alberta foreland basin (Fig. 5; Ross and Villeneuve, 2003; Lewis et al., 2010; Leier and Gehrels, 2011; Raines et al., 2013; Laskowski et al., 2013; Gehrels and Pecha, 2014).

A widely interpreted mechanism for westward transportation of Appalachian and Grenville detritus from eastern to western North America are continental-scale paleorivers (Rainbird et al., 1992; Dickinson and Gehrels, 2009; Dickinson et al., 2010; Benyon et al., 2014, 2016; Blum and Pecha, 2014; Finzel, 2014). In addition to Mesoproterozoic–Triassic passive-margin strata, a relevant sink for these sediments was vast eolian deposits in the southwestern United States, and these were likely an important source of detrital zircon grains to the Western Interior Basin (Dickinson and Gehrels, 2009; Leier and Gehrels, 2011; Laskowski et al., 2013; Raines et al., 2013).

Northern sources of detrital zircon grains exist in the North American passive margin (Lemieux et al., 2011; Hadlari et al., 2012; Gehrels and Pecha, 2014). However, these sources are discounted because the basin configuration was dominated by south to north sediment transport during the Mesozoic (Jackson, 1984; Leckie and Smith, 1992).

To aid provenance interpretations, five provenance groups were established with the aid of multidimensional scaling. This analysis relies on incorporating nearest-neighbor lines to define the groups of related spectra, and synthetic end-member populations to show important inputs to each group (Fig. 6).

Group 1: Trans-Hudson–Dominated Spectrum with Minor Cordilleran Arc

This group includes only the Monteith A member (Fig. 6). The detrital zircon spectrum is characterized by a prominent Trans-Hudson age at ca. 1850 Ma, with a subordinate mode at ca. 2080 Ma. This sample has only minor proportions of zircon grains from the other primary detrital zircon sources of North America (Fig. 4).

The Monteith A member is interpreted to have been derived principally from passive-margin strata via transversely oriented rivers that debouched directly from the incipient mountain belt (Kukulski et al., 2013a; Raines et al., 2013). Recycling of Neoproterozoic to Ordovician strata can account for the majority of zircon grains in the Monteith A member sample; the remainder could be derived from Devonian to Triassic strata, indicating that passive-margin strata were the primary sediment source for this unit (Fig. 5; Gehrels and Pecha, 2014). The youngest mode (200 Ma) indicates that paleorivers with catchments extending westward to the Cordilleran magmatic arc existed during deposition of the Monteith A member.

Group 2: Yavapai-Mazatzal/Trans-Hudson–Dominated Spectrum with Significant Cordilleran Arc

This provenance group includes samples from the Cadomin Formation, Bluesky Formation, and Dunvegan Formation (Fig. 6), which are all components of distinct tectonostratigraphic wedges (Table 1; Fig. 3). These samples are characterized by a dominant mode centered at ca. 1800 Ma. They all contain a significant proportion of zircon grains attributed to the Yavapai-Mazatzal Province (1620–1790 Ma), as well subordinate modes from every other major North American detrital zircon source. Group 2 is distinguished from group 1 by additional Cordilleran arc–derived zircon grains (Fig. 4).

This provenance group is interpreted to derive, in large part, from passive-margin strata of western North America. Like Group 1, these sediments are very similar to Neoproterozoic through Ordovician strata (Fig. 5; Gehrels and Pecha, 2014). The addition of Jurassic–Cretaceous zircon grains indicates increased input from the magmatic arc (Armstrong, 1988).

Group 3: Proterozoic-Dominated Spectrum

This provenance group includes the most samples and is the most internally diverse. It includes the Monteith C member (Raines et al., 2013), Gething Formation, Falher D Member, Cadotte Member, Cardium Formation, and Chinook Member (Fig. 6). These units are characterized by a more equal representation of the Proterozoic-age groups (i.e., Grenville, Granite-Rhyolite Province, Yavapai-Mazatzal, and Trans-Hudson; Fig. 4).

The samples in this group contain a wide spectrum of Proterozoic ages that are attributed to recycling from passive-margin strata (Gehrels and Pecha, 2014). However, unlike groups 1 and 2, the ca. 1800 Ma mode is subordinate. The proportions of different zircon populations are more similar to Devonian–Triassic strata of the passive margin than to Neoproterozoic–Ordovician deposits (Fig. 5; Gehrels and Pecha, 2014). The presence of Mesozoic ages indicates incorporation of the magmatic arc in the catchments of each unit.

Group 4: Paleozoic–Neoproterozoic–Dominated Spectrum

This provenance group includes the Monteith C member, Monteith B member, and Falher A Member and is characterized by the relative importance of Paleozoic and Neoproterozoic detrital zircon grains that are interpreted to be associated with the Appalachian orogen (285–850 Ma; Fig. 4). This group is also characterized by relative importance of detrital zircons attributed to the Grenville orogen over detrital zircons from the Trans-Hudson orogen (Fig. 6B). The Cretaceous Falher A Member has a higher proportion of Cordilleran zircon grains than the Jurassic Monteith Formation samples in this group, indicating increased availability of Cordilleran arc detritus.

Similar Paleozoic and Neoproterozoic detrital zircon ages occur in Triassic passive-margin strata of western North America and in eolianites in the southwestern United States (Fig. 5; Dickinson and Gehrels, 2009; Leier and Gehrels, 2011; Gehrels and Pecha, 2014; Golding et al., 2016). On the basis of a more complete data set including paleoflow measurements and petrography, the Monteith C and Monteith B members were interpreted to have received significant contributions from basin-axial paleorivers with catchments including the southwestern United States and/or the Appalachian system directly (Fig. 5; Raines et al., 2013). An axial river source is the favored interpretation for group 4 based on various investigations (e.g., Hamblin and Walker, 1979). This axial paleoriver system may have reactivated during deposition of the Falher A Member, which is consistent with paleogeographic interpretations (Leckie and Smith, 1992).

Group 5: Cordilleran Arc–Dominated Spectrum

This provenance group includes only the Notikewin Member, where 75% of the ages in this sample were derived from the Cordillera (<252 Ma; Fig. 4). This is 52% more than the next nearest sample (Cardium Formation).

The Notikewin Member is interpreted to have derived from paleorivers with catchments that included the Cordilleran magmatic arc. Jurassic detrital zircon grains (ca. 165 Ma) arrived in the basin as early as the Late Jurassic (Monteith C member), but they did not dominate until deposition of the Notikewin Member in the Albian. Therefore, deposition of this member is likely associated with a time of major Middle Jurassic pluton unroofing. Plutons of this age are common in the Omineca belt of interior British Columbia (Archibald et al., 1983; Armstrong, 1988).

DeCelles et al. (2009) hypothesized linkages between upper- and lower-plate processes. In their model, underthrusted retroarc crust melts and forms a dense eclogitic root that founders into the mantle. Following this foundering, a period of increased magmatism occurs, termed a high-flux episode (Fig. 7). This increase in magmatism is associated with plateau uplift and extension in the interior of the orogen, changing the taper angle of the orogenic wedge and driving a period of shortening in the retroarc thrust belt as the plateau collapses and spreads laterally. A period of quiescence occurs when the taper of the orogenic wedge becomes subcritical; subsequently, the crust is pulled down during the formation of the next eclogitic root (DeCelles et al., 2009). These processes should be manifested in the fill of a foreland basin in four ways: (1) increased sediment supply during periods of uplift; (2) propagation of a flexural wave creating basin subsidence during thrust loading; (3) evolving provenance from migrating drainage divides driven by cyclical uplift; and (4) evolving provenance linked to evolving sediment routing systems (DeCelles et al., 2009). The formation of unconformities is not explicitly predicted in this model; however, we propose that periods between high-flux episodes may be linked with times of uplift in the basin in response to denudation of thrust loads.

We compiled detrital zircon ages younger than 200 Ma from foreland basin strata of Alberta, Montana, and British Columbia (Fig. 8). A key interpretation is that prominent detrital zircon age modes for the basin correspond to major episodes of magmatism in the adjacent arc, which correspond to high-flux episodes (DeCelles et al., 2009; Laskowski et al., 2013). There are three major modes in the probability density function at ca. 165, 115, and 74 Ma (Fig. 8). The proportion of zircon grains in each mode is interpreted as a bias toward sampling of older strata that predate the younger modes, and the availability of the older detrital zircon grains to young and old strata; therefore, the prominence of the modes does not represent intensity of an individual magmatic episode. Other factors that may contribute to bias are differences in zircon fertility (e.g., Dickinson, 2008), depth of magma emplacement, and preservation bias (Hawkesworth et al., 2009; Spencer et al., 2015b). Selective preservation of zircon from the collisional phase of supercontinent cycles and preferential destruction from subduction and rifting phases have been described by Hawkesworth et al. (2009) and Spencer et al. (2015b). While preservation bias could potentially explain episodic peaks in the probability density plots of Cordilleran detrital zircon peaks from the foreland basin, the recurrence interval of these peaks (<50 m.y.) is not on the same scale as the supercontinent cycles investigated by previous workers.

Compiling detrital zircon ages from the basin is complementary to bedrock mapping for understanding the magmatic history of the Canadian Cordillera. The sedimentary system is arguably a more efficient way of sampling the arc because river systems act as a random sampling mechanism and concentrate a large number of zircon grains that can be analyzed to produce a robust data set. Using the basin as a means to sample the arc is also advantageous because detrital zircon grains in the basin were derived from plutons and volcanic deposits that may have been completely eroded.

Ar-Ar fault gouge cooling ages from the fold-and-thrust belt and K-Ar dates from plutons in the Kootenay arc of south-central British Columbia are included with the detrital zircon data and show a general temporal association with the episodic modes in the detrital zircon data (Fig. 8; Archibald et al., 1983; Pana and van der Pluijm, 2015). This association indicates that the proposed high-flux episodes were periods of orogen denudation and thrusting in the retroarc. The durations of major unconformities are also plotted to highlight periods of uplift in the basin. With this framework of magmatic pulses, the provenance evolution of the stratigraphy at Grande Cache and the subdivision of the stratigraphy into tectonostratigraphic packages are given broader context.

The first major high-flux episode is proposed to have occurred in the Canadian Cordillera at ca. 165 Ma. According to fault gouge cooling ages, this high-flux episode was followed by a period of thrusting as predicted by the model (Fig. 8; DeCelles et al., 2009; Pana and van der Pluijm, 2015). We interpret that uplift following the foundering of the crustal root led to increased sediment supply and the infilling of the basin as recorded by the first tectonostratigraphic wedge (Fig. 1). The provenance during this time is interpreted to have evolved from basin-axial input rich in Appalachian detrital zircon grains recycled from the southwestern United States, culminating in provenance dominated by recycled zircon grains derived from passive-margin rocks to the west of the study area that were uplifted during this first episode of deformation (Figs. 4 and 6). This change reflects evolution from an underfilled marine basin with low sediment input from the adjacent orogen to a basin dominated by westerly derived sediment (Miles et al., 2012). The transition from axial river systems to transverse river systems in a foreland basin has been attributed to decreased tectonic loading and uplift in the orogen caused by erosion-driven isostatic uplift (Burbank, 1992). This is consistent with the interpreted provenance evolution of the Montieth Formation (Raines et al., 2013).

From 145 to 125 Ma, the sub-Cretaceous unconformity formed, attributed to basin uplift as orogenic activity waned (Heller et al., 1988; Cant and Stockmal, 1989; Ross et al., 2005). The sub-Cretaceous unconformity corresponds to a time of decreased magmatism inferred from the rarity of zircon grains of this age in the foreland basin (Fig. 8). This period of basin uplift transitioned into the deposition of the Bullhead Group, which does not thicken abruptly westward and is considered by some to record a period of erosion and reworking, as opposed to a phase of high subsidence and elevated sedimentation (Cant and Stockmal, 1989; Heller and Paola, 1989; Hayes et al., 1994). The provenance of the Cadomin and Gething Formations is interpreted to be recycled passive-margin deposits, though the basin-scale sediment routing system was characterized by axial river systems during this time (Jackson, 1984; Hayes et al., 1994).

The second high-flux episode is proposed to have occurred at ca. 115 Ma (Fig. 8). This is also associated with fault gouge dates overlapping the waning phase of apparent magmatic activity (Pana and van der Pluijm, 2015). Following this event, a major period of basin subsidence and sedimentation occurred, associated with deposition of the Spirit River Formation (Leckie and Smith, 1992). Stott (1984) estimated that sedimentation rates increased during deposition of the Fort St. John Group, based on biostratigraphic control and measured stratigraphic thicknesses. The provenance during this interval evolved from a North American passive-margin–derived signature (i.e., Falher D Member), to one that likely represents axial input (i.e., Falher A Member), to a Mesozoic Cordilleran-dominated signature (i.e., Notikewin Member; Fig. 4). The preponderance of Cordilleran-derived Mesozoic zircon grains in the Notikewin Member is attributed to denudation in the orogen sufficient to unroof ca. 165 Ma plutons emplaced in the Omineca belt during the first high-flux episode. The drainage divide is interpreted to have migrated far to the west during deposition of the Notikewin Member.

The ca. 115 Ma mode is associated with a subordinate mode at ca. 104 Ma potentially indicating that Aptian–Albian magmatic activity extended over a protracted time period (Fig. 8). The subordinate mode overlaps with the duration of the sub-Paddy unconformity (Stelck and Leckie, 1990). This makes the association between magmatic lulls and unconformities less evident than in the case of the sub-Cretaceous unconformity.

The final phases of sedimentation in the study area are more difficult to link to the magmatic flux proxy. The age of the last high-flux episode at ca. 74 Ma postdates the strata of this study. The Cadotte-Chinook interval could either be related to a continuation of the second tectonostratigraphic wedge or a time of tectonic quiescence (cf. Cant and Stockmal, 1989; Pana and van der Pluijm, 2015). The detrital zircon provenance of the Cadotte-Chinook interval shows an evolution from broad, passive-margin spectra to those increasingly dominated by Cordilleran detritus. This pattern is similar to that of the second tectonostratigraphic wedge that culminated in deposition of the Notikewin Member (Fig. 4). The cumulative stratigraphic thickness of the Codotte-Chinook interval is ∼1200 m, indicating significant subsidence (Fig. 3). In general, these units are characterized by shale-dominated sedimentation in the basin punctuated by progradational cycles of coarser deltaic and shoreface sediments (Lerand, 1983; Smith et al., 1984; Walker, 1983; Rahmani and Smith, 1988; Leckie, 1989; Leckie and Smith, 1992; Beaumont et al., 1993; Bhattacharya, 1994; Plint, 2000; Collom, 2001; Shank and Plint, 2014; Benham and Collom, 2012). Sediment supply estimates are low for the Dunvegan, Cardium, and Chinook units in comparison to the Fort St. John Group and Monteith Formation–Upper Minnes Group (Stott, 1984). A possible explanation of these data is that periods between high-flux episodes are predicted by the model to be characterized by mountain belts with less elevation, possibly contributing to reduced sedimentation in the adjacent basin (DeCelles et al., 2009).

The Campanian magmatic event (ca. 74 Ma) is more closely associated in time with Brazeau-Paskapoo sedimentation of the Saunders Group (Figs. 1 and 3). This interval is dominated by fluvial sediments and elevated sedimentation rates, although it does not outcrop in the vicinity of Grande Cache and is therefore not included in this study (Stott, 1984; Leckie and Smith, 1992). The Saunders Group may represent the final phase of high sediment input and basin filling following the proposed Late Cretaceous high-flux episode. This package of sediment has been described as the thickest accumulation of sediment in the basin deposited during the time of maximum thrust shortening (Ross et al., 2005).

The analysis presented here understates the potential significance of terrane accretion events in the Canadian Cordillera. Both terrane accretion events and high-flux episodes occurred episodically, and each analysis has shortcomings when comparing the timing of Cordilleran events to the stratigraphic record (Figs. 1 and 8; Price et al., 1981; Cant and Stockmal, 1989).

The cyclical model is attractive because it incorporates a more diverse suite of processes that are interpreted to operate in Cordilleran systems (DeCelles et al., 2009). The cyclical model is also advantageous because the absolute ages of magmatic events are more precisely constrained than the timing of accretion events, allowing for a more rigorous comparison of orogenic events and the preserved foreland stratigraphy (Price et al., 1981; Cant and Stockmal, 1989).

The terrane accretion model allows for a more straightforward understanding of load emplacement and flexural response of the lithosphere (Stockmal and Beaumont, 1987; Cant and Stockmal, 1989). Critiques of the terrane accretion model include: (1) mismatched scales between accreted terranes and sediment volumes in the foreland; (2) the easterly transmission of lithospheric stresses associated with accretion events that occurred far outboard of the foreland basin (i.e., the Bridge River terrane, ∼500 km; Coney et al., 1980); and (3) the possibility of large strike-slip reorganization of the orogenic collage (Cowan et al., 1997). These issues necessitated the caveat that small and far outboard accretion events do not influence basin subsidence directly but act as drivers of deformation toward the retroarc foreland basin, therefore acting as a mechanism for deposition of tectonostratigraphic wedges (Stockmal and Beaumont, 1987; Cant and Stockmal, 1989; Beaumont et al., 1993). This caveat may be validated by observations in the northern North American Cordillera, where the Yakutat block is accreting to the Alaska margin and driving shortening over 800 km away in the backarc (Mazzotti and Hyndman, 2002). There are no data to suggest that this collision is contributing significant sediment to the foreland basin.

More detailed studies are needed to assess the relationship between terrane accretion events and other orogenic processes. The high-flux episodes occurred during and in-between terrane accretion events, and there are more terrane accretion events reported for the Canadian Cordillera than there are high-flux episodes (Figs. 1 and 8). Therefore, the argument that they are linked is problematic. This may necessitate favoring one model over the other, or constructing a hybrid model that integrates cyclical uplift and episodic terrane accretion to more completely honor the stratigraphic record.

The stratigraphy of the Alberta foreland basin has been subdivided into tectonostratigraphic wedges, which provide the context for analysis of provenance evolution in Late Jurassic–Late Cretaceous foreland basin fill that crops out in the central Alberta Foothills. Each significant sandstone unit in the basin from the Late Jurassic to the middle of the Late Cretaceous contains zircon spectra that span the magmatic history of North America. The zircon spectra can be separated into five provenance groups using multidimensional scaling. Groups 1–3 are dominated by recycled sediment from the passive margin of western North America. Group 4 consists of sediment recycled from southern sources delivered via basin-axial fluvial transport. Group 5 is dominated by Jurassic zircon grains from the Cordilleran magmatic arc. These groups do not show a correlation to lithostratigraphic position or previously determined tectonostratigraphic wedge.

Mesozoic-age detrital zircon grains from Alberta, British Columbia, and Montana show three major modes at ca. 165, 115 and 74 Ma. A model links these modes to cyclical high-flux episodes of magmatism triggered by foundering of an eclogitic root in the orogenic belt. The first two high-flux episodes are linked to deposition of tectonostratigraphic wedges 1 and 2 of the Minnes Group and the Fort St. John Group, respectively. These sediments are separated by an intervening unconformity (10–20 m.y.) and magmatic lull that are interpreted as a time of basin uplift. The uppermost stratigraphy in this study records a package of sediment that predates the last high-flux episode. Estimates of sedimentation rates suggest that these units were deposited during a period of high subsidence and low sediment supply. The final high-flux episode is linked to a particularly thick stratigraphic succession of nonmarine deposits, which was not an emphasis in this study.

This analysis shows the feasibility of linking underplating, eclogitic root foundering, episodic magmatism and uplift, and collapse of the orogen with foreland basin fill, informing a more integrated understanding of cyclicity in Cordilleran orogenic systems. The cyclical model presented here for the organization of the foreland basin stratigraphy deemphasizes the effect of terrane accretion as the key driver of Alberta foreland basin deposition, although it is perhaps likely that the strata record a complex interplay between cyclical orogen development and sequential terrane accretion.

Funding for this research was generously provided by Shell Canada (undergraduate thesis of R. van Drecht) and a Natural Sciences and Engineering Research Council Discovery grant to S. Hubbard. Dale Leckie provided helpful discussions. Reviews by Christopher Spencer, Andrew Miall, and Science Editor Kurt Stüwe improved the clarity of the manuscript and are greatly appreciated.

1GSA Data Repository Item 2016148, a detailed description of the analytical methods used and U-Pb isotopic data and ages of detrital zircon grains, is available at www.geosociety.org/pubs/ft2016.htm, or on request from editing@geosociety.org.