sexta-feira, 23 de fevereiro de 2018

Fossil and genomic evidence constrains the timing of bison arrival in North America

Duane Froese, Mathias Stiller, Peter D. Heintzman, Alberto V. Reyes, Grant D. Zazula, André E. R. Soares, Matthias Meyer, Elizabeth Hall, Britta J. L. Jensen, Lee J. Arnold, Ross D. E. MacPhee and Beth Shapiro
  1. Edited by Donald K. Grayson, University of Washington, Seattle, WA, and approved February 3, 2017 (received for review December 20, 2016)

Significance

The appearance of bison in North America is both ecologically and paleontologically significant. We analyzed mitochondrial DNA from the oldest known North American bison fossils to reveal that bison were present in northern North America by 195–135 thousand y ago, having entered from Asia via the Bering Land Bridge. After their arrival, bison quickly colonized much of the rest of the continent, where they rapidly diversified phenotypically, producing, for example, the giant long-horned morphotype Bison latifrons during the last interglaciation.

Abstract

The arrival of bison in North America marks one of the most successful large-mammal dispersals from Asia within the last million years, yet the timing and nature of this event remain poorly determined. Here, we used a combined paleontological and paleogenomic approach to provide a robust timeline for the entry and subsequent evolution of bison within North America. We characterized two fossil-rich localities in Canada’s Yukon and identified the oldest well-constrained bison fossil in North America, a 130,000-y-old steppe bison, Bison cf. priscus. We extracted and sequenced mitochondrial genomes from both this bison and from the remains of a recently discovered, ∼120,000-y-old giant long-horned bison, Bison latifrons, from Snowmass, Colorado. We analyzed these and 44 other bison mitogenomes with ages that span the Late Pleistocene, and identified two waves of bison dispersal into North America from Asia, the earliest of which occurred ∼195–135 thousand y ago and preceded the morphological diversification of North American bison, and the second of which occurred during the Late Pleistocene, ∼45–21 thousand y ago. This chronological arc establishes that bison first entered North America during the sea level lowstand accompanying marine isotope stage 6, rejecting earlier records of bison in North America. After their invasion, bison rapidly colonized North America during the last interglaciation, spreading from Alaska through continental North America; they have been continuously resident since then.
The invasion of bison (Bison) from Asia across the Bering Isthmus profoundly affected the North American faunal community. Bison, or American buffalo, are large-bodied, aggressive, and highly fecund. Following their establishment in North America, bison rapidly became the most important competitor for forage within the established large mammal community (1). Early North American bison were morphologically and ecologically diverse (2, 3). In addition to extant Bison bison, taxa traditionally recognized in systematic treatments include the steppe bison (Bison priscus), which first colonized northwestern North America from Asia, and the giant long-horned bison (Bison latifrons) of the central and southern continent (Figs. 1A and 2). The latter is the largest bison known. It inhabited woodlands and forest openings through much of the continental United States and southern Canada; however, their fossils have not been found in northern Canada or Alaska. Bison eventually became an important hunting resource for Indigenous North Americans (4) and remain an icon of the American plains (5).
Fig. 2.
Reconstructions of bison skulls based on fossils attributed to (A) a giant long-horned bison, B. latifrons; (B) a Late Pleistocene steppe bison, B. priscus; and (C) a present-day B. bison. Giant long-horned bison were significantly larger than present-day bison; adult males may have weighed in excess of 2,000 kg, which is twice as large as present-day bison, and had horns that spanned as much as 2.2 m (57, 58).
At the time of bison arrival, the North American megaherbivore grazing community was dominated by mammoths (Mammuthus) and caballine equids (Equus). Equids have a deep evolutionary history in North America that is closely associated with the rise of grasslands during the Early Miocene (∼18 Ma) (6). Mammoths dispersed from Asia to North America during the Early Pleistocene (∼1.35 Ma), becoming the continent’s largest-bodied obligatory grazer (7). The subsequent arrival of Bison markedly affected a faunal community dominated by Equus and Mammuthus, but when that process actually began has been difficult to determine.
Paleontologically, Bison is the index taxon for the Rancholabrean, the final North American Land Mammal Age (810). Land Mammal ages are important because, in the absence of other chronological data, they provide a means to infer the age of a locality based on taxonomic assemblages. This assumes, however, that the first appearance datum of the index taxon can be reliably tied to a specific interval, which is surprisingly not the case for bison (8). Here, we used a combined paleontological and paleogenomic approach to establish the timing of bison entry into North America.
Models of the timing of Bison entry into North America range, based on fossil occurrences, from arrival during the Late Pliocene or Early Pleistocene, approximately 2–3 Ma, to the Late Pleistocene (8). The oldest of these dates are based on fossil sites in Florida (11) and Alaska (12) that have now been shown to be of poor stratigraphic and chronologic integrity (8, 13, 14). A bison astragalus dating to 290–230 thousand y before present (kyBP), or marine isotope stage (MIS) 7, was recovered from South Carolina (15), but the reliability of this age and its association with the fossils have also been questioned (16, 17). Bison fossils are absent in well-dated Middle Pleistocene (780–130 kyBP) localities from both central and northern North America (8). In Kansas, bison remains are all younger than the Lava Creek B Ash, which has a maximum constraining age of 640 kyBP (18). Bison fossils are also absent from the rich Sheridanian fauna “Equus beds” of Nebraska (19), which are overlain by the Loveland loess, dated regionally using optically stimulated luminescence to 180–130 kyBP, or MIS 6 (20). In Yukon, a fossil assemblage in stratigraphic association with the Middle Pleistocene Gold Run tephra (735 ± 88 kyBP) includes horse (Equus), proboscideans (Mammuthus), sheep (Ovis), and biostratigraphically diagnostic microtine rodent species, but bison fossils are absent (2123). The oldest site in the midcontinental United States possessing both bison and a firm chronology is near Snowmass, Colorado, where bison fossils were recovered within sediments associated with the last interglaciation (MIS 5d), thus dating to ∼120 kyBP (24). Bison fossils also occur at American Falls, Idaho, associated with a lava-dammed lake that dates to ∼72 kyBP (8, 25). Bison from both of these sites represent the giant long-horned bison, B. latifrons (24, 25), generally considered to be the earliest form of bison developed in the continental United States.
To establish a reliable first appearance datum for bison in North America, we first characterized the in situ fossil assemblage and chronology of two, well-dated fossil-rich localities in the Old Crow area of northern Yukon, Canada: CRH 11 and Ch’ijee’s Bluff (Figs. 1 B and C and 3). Next, we isolated and sequenced mitochondrial genomes from two of the oldest Bison fossils yet identified: a partial metacarpal found at Ch’ijee’s Bluff that dates to ∼130 kyBP (Fig. 3C), and a humerus from the site near Snowmass recovered within a layer dated to ∼120 kyBP (Fig. 1D). These two bison were identified, based on their size and geographic location, as a steppe bison and a giant long-horned bison (24), respectively. We then used a coalescent-based approach to infer a mitochondrial genealogy for these and 44 other Late Pleistocene and Holocene bison, taking advantage of an approach to calibrate a molecular clock within a coalescent framework using the ages of each sampled bison (26). This approach allows an estimate of the age of relevant nodes in the bison mitochondrial genealogy, making it possible to test hypotheses not only about the timing of bison entry into North America but also about the relationship between these morphologically distinct bison forms.
Fig. 3.
Features of the Ch’ijee’s Bluff locality. (A) the Old Crow tephra (124 ± 10 kyBP; UA1206) highlighted by the white arrow, (B) the stratigraphic setting of the Old Crow tephra, bison metacarpal YG 264.1, and the MIS 5e forest bed, and (C) the in situ metacarpal was found several centimeters beneath the prominent MIS 5e forest bed and ∼125 cm above Old Crow tephra (the latter is not shown). The stratigraphy indicates a latest MIS 6 age for YG 264.1.

Results and Discussion

Characterizing Two Fossil Assemblages in Yukon, Canada.

To assess the chronology of bison presence in high-latitude northwest North America, we first characterized the in situ fossil assemblage and chronology at CRH 11 (27, 28) and Ch’ijee’s Bluff (Fig. 3). The chronologies of both of these sites rely heavily on identification of volcanic ash layers, or tephra, in sediment exposures. Individual tephra are readily characterized geochemically and, when combined with other correlative indicators—such as stratigraphy or paleoecology—can provide isochronous stratigraphic markers across a region (29). The Old Crow tephra has an isothermal plateau glass fission-track age of 124 ± 10 kyBP at the 1σ confidence level (29). This age determination spans the transition from the late MIS 6 glacial to the MIS 5e interglaciation, but paleoecological evidence for cool climate conditions during tephra deposition indicates a late MIS 6 age (29, 30). The Old Crow tephra thus provides a useful marker horizon for this study: interglacial deposits below the tephra must date to MIS 7 or older, whereas sediment above the tephra but below the prominent interglacial deposits represents a narrow time interval between latest MIS 6 and the beginning of MIS 5e (∼135 to ∼125 kyBP).
CRH 11 is one of the classic localities for Quaternary paleontology in North American Beringia. This bluff, on the left bank of the Old Crow River (67° 49′ N, 139° 51′ W) comprises ∼30 m of silt and sand that is locally organic-rich (28). We recovered 294 vertebrate fossils from the lowest bone-bearing unit at CRH 11, including specimens of woolly mammoth (Mammuthus primigenius), horse (Equus sp.), caribou (Rangifer tarandus), giant beaver (Castoroides ohioensis), beaver (Castor canadensis), wolverine (Gulo gulo), and Jefferson’s ground sloth (Megalonyx jeffersonii) (for a complete list, SI Text). Bison fossils are absent, consistent with earlier in situ assemblages recovered from the site (27). The Old Crow tephra is present 14 m above the in situ fossils (Fig. S1), establishing that the fossiliferous sediments must be older than the last interglaciation (pre-MIS 5e) and, based on paleoecology, date from the penultimate interglaciation, MIS 7 (28). To confirm this age, we performed direct single-grain optically stimulated luminescence dating of these sediments (SI Text), which gave a weighted mean age of 208 ± 6 kyBP (Figs. S2S4), consistent with a MIS 7 age. Bison fossils were absent not only from this stratigraphic level at CRH 11, but also from nearby sites of comparable age (27), suggesting that bison were not present in Yukon before MIS 6.
Fig. S1.
Major element geochemistry of Old Crow tephra at CRH 11 and Ch’ijee’s Bluff.
Fig. S2.
OSL dose-recovery test results for sample CRH 11-1. (A) Ratios of recovered-to-given dose versus PH1 temperature (held for 10 s) for ∼400-grain aliquots. The natural OSL signals of the multigrain aliquots were optically bleached with two 1-ks blue LED illuminations at ambient temperature, each separated by a 10-ks pause. A known dose of 100 Gy was then administered to each aliquot and the SAR procedure detailed in table S2 of Arnold et al. (41) was subsequently used to estimate this dose. A fixed PH2 of 160 °C for 10 s was applied in the SAR procedure. (B) Radial plot showing ratios of recovered-to-given dose obtained for individual quartz grains using a PH1 of 160 °C for 10 s and a PH2 of 160 °C for 10 s. The single-grain natural OSL signals were bleached using the same procedure outlined above and the administered doses were subsequently recovered using the single-grain SAR procedure shown in Table S1. The gray shaded region on the radial plot is centered on the administered dose for each grain (sample average = 315 Gy; although this amount varied from 265 to 353 Gy for individual grains, because of spatial variations in the dose rate of the β source). Individual De values that fall within the shaded region are consistent with the administered dose at 2σ.
Fig. S3.
Representative OSL decay and dose–response curves for individual quartz grains from the CRH 11 samples. Sensitivity-corrected dose–response curves were constructed using the first 0.17 s of each green laser stimulation after subtracting a mean background count obtained from the last 0.25 s of the OSL signal. (A) Quartz grain from CRH 11-2 with a relatively bright OSL signal (Tn >500 counts/0.17 s). (B) Quartz grain from sample CRH 11-3 with an average OSL signal (Tn = ∼200–300 counts/0.17 s). (C) Quartz grain from CRH 11-2 with a particularly high characteristic saturation dose (D0) value of >350 Gy. In the inset plots, the open circles on the y axis denote the sensitivity-corrected natural OSL signals, and the sensitivity-corrected regenerated OSL signals are shown as filled circles.
Fig. S4.
Single-grain De distributions for the four CRH 11 OSL samples, shown as radial plots (A–D). The De value for each grain is read by drawing a line from the origin of the y axis (Standardized Estimate), through the data point of interest, to the radial axis (plotted on a log scale) on the right-hand side; the point of intersection is the De (in Gy). The measurement error on this De is obtained by extending a line vertically to the x axis, where the point of intersection is the relative SE (shown as a percentage of the De value) and its reciprocal (Precision). In radial plots, the most precise estimates fall furthest to the right, and the least precise estimates fall furthest to the left. Here, the gray bands are centered on the weighted mean De values used to calculate the OSL ages, which we estimated using the CAM of Galbraith et al. (43).
Bison fossils are present, however, at nearby Ch’ijee’s Bluff 40 km southwest (67° 29′ N, 139° 56′ W), where 50 m of unconsolidated Late Cenozoic sediment are exposed along a prominent 4-km-long cut bank on the Porcupine River (Figs. 1C and 3). The Old Crow tephra (29) is present across much of Ch’ijee’s Bluff (Fig. 3A), where it underlies a prominent bed of dark-brown macrofossil-rich organic silt. Paleoecological indicators within this organic layer represent a closed boreal forest and a warmer than present climate, indicating a MIS 5e age (30, 31). In support of this age and the stability of the organic layer, logs and organic detritus from the organic bed all have nonfinite 14C ages (30, 31). Although most bison fossils at the site are detrital, we recovered a single in situ bison metacarpal (YG 264.1) 125 cm above the Old Crow tephra and directly below the prominent MIS 5e woody peat bed (Fig. 3 B and C). The tephra- and bone-bearing silt containing the fossil is sharply overlain by the organic-rich bed, an unconformity related to the thawing of permafrost during the last interglaciation that is common across northwestern Canada and Alaska (30). Given this stratigraphic context, we conclude that the bison fossil dates to ∼130 kyBP: it is younger than late MIS 6 (the age of Old Crow tephra), and older than the onset of the MIS 5e interglaciation. This is the oldest reliably dated fossil evidence of bison in North America (cf. ref. 8).

The Oldest Fossil Bison in North America.

Bison fossils are common throughout Yukon, Alaska, and Siberia. All are medium-horned bison, most commonly referred to as steppe bison, B. priscus (2), and published genetic data are consistent with this interpretation (3). The discovery of giant long-horned bison at a last interglacial site (MIS 5 sensu lato) near Snowmass, Colorado, establishes the presence of a morphologically distinct bison in continental North America (Fig. 2A). These giant long-horned forms have never been recovered from northern locales between Yukon and Siberia (2). The age of the Snowmass site places the long-horned bison fossil slightly younger than the bison from Ch’ijee’s Bluff.
Because of their antiquity and the poor preservation conditions of continental compared with northern (permafrost) localities, giant long-horned bison fossils have thus far failed to yield usable DNA. However, recent advances in paleogenomics have expanded the range of fossils from which DNA can be recovered (32). Capitalizing on these, we used a hybridization capture approach to enrich for bison mitochondrial DNA from both the Ch’ijee’s Bluff bison and from a giant long-horned bison from Snowmass (DMNH EPV.67609) that dates to the last interglaciation, ∼120 kyBP (24). We recovered a complete mitochondrial genome (159× coverage) from the Ch’ijee’s Bluff bison and a near-complete mitochondrial genome (6.6× mean coverage) from the Snowmass bison (Materials and Methods). We assembled complete mitochondrial genomes from an additional 6 Siberian and 26 North American bison ranging in age from ∼0.4–45 kyBP (Dataset S1). We then estimated the evolutionary relationship between these and 10 present-day American (33, 34) and an ancient Siberian bison (35), using stratigraphic ages and radiocarbon dates to inform the molecular clock (Materials and Methods and SI Text).
Both the Ch’ijee’s Bluff bison and the Snowmass bison mitochondrial lineages fall near the root of sampled bison mitochondrial diversity, indicating that both bison were early descendants of the first bison dispersing into North America (Fig. 1E). North American bison share a common maternal ancestor 195–135 kyBP (Fig. 1E, node 1, SI Text, and Fig. S5), consistent with the MIS 6 glaciation (Fig. 1E), and precluding models for a significantly older bison presence in North America. This timing is coincident with an interval of reduced eustatic sea level that would have enabled interchanges across the Bering Isthmus (36, 37). We also identify a second, later dispersal of bison from Asia into North America during the Late Pleistocene, ∼45–21 kyBP (Fig. 1E, node 2, SI Text, and Fig. S5), within a period of lowered sea level during the last glaciation (36).
Fig. S5.
Maximum clade credibility trees resulting from the Bayesian analysis of the reduced (A) and full (B) mitochondrial genome data sets. The reduced dataset (A) is also depicted in Fig. 1E. Colors are as in Fig. 1. Numbers along the branches are Bayesian posterior probability scores for each clade. Bars represent 95% highest posterior probability density intervals for node heights and are reported for nodes with a posterior probability score of >0.95 (A) or >0.85 (B). Tip labels follow the convention of sampleID_locality/species_age as in Dataset S1.
The Rancholabrean is the most recent of the North American Land Mammal Ages, and has long been defined by the presence of bison in continental records (9). However, Bell et al. (8) argued that this North American Land Mammal Age should only apply to localities or faunules recovered from south of 55° North. Their reasoning was that, given its proximity to the Bering Land Bridge and eastern Eurasia, the northern part of the continent required a distinct chronology because of the potential for faunal mixing. Our results, however, show temporal and genetic affinity between the arrival of bison in northwestern Canada and their dispersal further south. The close genetic relationship between maternal lineages found in the earliest northern bison and the earliest continental bison argues for a rapid expansion of bison across the continent in a period of approximately 20,000 y between late MIS 6 and MIS 5d. These records also demonstrate the rapid phenotypic change from northern forms of bison (e.g., B. priscus, B. alaskensis) found in Siberia through Yukon, to B. latifrons in the continental United States.
The integration of independent geochronological data with faunal collections and a molecular dating approach constrains the history and dynamics of bison dispersal into North America. These complementary approaches provide a remarkably consistent picture of this grazer as it entered the continent during the sea level lowstand accompanying MIS 6 and spread from Alaska through the continental United States. The rapid dispersal and success of bison in North America make a strong case for bison as an index taxon for the Rancholabrean at a continental scale. Although full nuclear genomic resources for bison are not yet available, these well preserved specimens will be important to future work to better understand the genetic basis for the remarkable phenotypic variability in early North American bison. Given their relatively shallow history and success in North American ecosystems, the entry of bison stands with human arrival as one of the most successful mammalian dispersals into North America during the last million years.

SI Text

Stratigraphic Information and Chronological Contexts of CRH 11 and Ch’ijee’s Bluff.

Old Crow tephra.

Chronology at Ch’ijee’s Bluff and CRH 11 relies on identification of tephra in sediment exposures. Individual tephra are readily characterized geochemically and, when combined with other correlative indicators, such as stratigraphy or paleoecology, can provide isochronous stratigraphic markers across a region.
Old Crow tephra, a rhyolite, is one of the most prominent of the numerous tephra preserved in Late Cenozoic sediments of eastern Beringia. Old Crow tephra is found throughout Alaska and western Yukon, with an estimated eruptive volume of ∼200 km3 (29). Old Crow tephra has an isothermal plateau glass fission-track age of 124 ± 10 kyBP (weighted-mean of four age determinations at 1σ) (29), which is consistent with independent thermoluminescence (59) and infrared-stimulated luminescence (60) ages on bracketing loess near Fairbanks, Alaska. Throughout its range of distribution, Old Crow tephra is present in relatively organic-poor, massive to crudely stratified silt, with paleoecological indicators for a cool climate at the time of deposition (29, 38). The tephra is up to 5 m below woody or peaty organic-rich silt, with pollen, plant, and insect macrofossil indicators for boreal forest conditions and climate as warm or warmer than present that is ascribed to the peak of the last interglaciation, MIS 5e (31, 6164). These age determinations, together with regional stratigraphic and paleoecological data, are most indicatative that the Old Crow tephra was deposited during late MIS 6 (2931).
Old Crow tephra provides a robust chronostratigraphic marker horizon for this study: interglacial deposits below the tephra must predate MIS 5, whereas sediment above the tephra but below the prominent last interglacial wood-rich/peaty silt represents a relatively narrow time window during latest MIS 6. Because the bison fossil was recovered from this narrow context, we ascribe a latest MIS 6 age to it of ∼130 kyBP. This age falls within the 1σ error range of the Old Crow tephra of 124 ± 10 kyBP (29, 30).
Field identification of Old Crow tephra was confirmed by electron probe microanalysis of glass shard major element geochemistry (Fig. S1) (28, 30). Briefly, sieved glass separates were mounted in acrylic pucks and analyzed on a JEOL 8900 microprobe at University of Alberta by wavelength dispersive spectrometry. Analytical conditions were a beam diameter of 10 μm, 15-keV accelerating voltage, and 6-nA beam current; secondary standards (Lipari obsidian and a reference sample of Old Crow tephra) were used for quality control and to minimize potential variation caused by differences in calibration and standardization between analytical sessions.

CRH 11.

Stratigraphy.
CRH 11, on the left bank of the Old Crow River (67° 49′ N, 139° 51′ W), is one of the classic localities for Quaternary paleontology in eastern Beringia (27, 65, 66). The bluff comprises ∼30 m of silt and sand, locally organic-rich, as described over 3 d of fieldwork in 2008 (28). The lowermost ∼2 m of massive sandy-silt are unconformably overlain by ∼15 m of low-angle cross-bedded sand and local gravel, grading upward to laminated silt and sand with abundant plant litter, and finally massive sandy silt at ∼17 m above river level. This unit is overlain by 3 m of massive sandy silt with ice-wedge casts and disseminated pods of Old Crow tephra (28). The tephra-bearing unit is overlain by ∼7 m of stratified sand with periglacial structures and organic interbeds, capped by a 1-m-thick peat. The uppermost 3 m of sediment at CRH 11 are stratified sand and silt with organic detritus, grading to dark gray clays associated with Glacial Lake Old Crow (6769).
Chronology.
Chronology for the CRH 11 site is based on the presence of the Old Crow tephra at 18 m above river level (28) and OSL dating of sands from the lower bone-bearing unit up to 3 m above river level Because the Old Crow tephra marks latest MIS 6 time, the interglacial unit underlying the Old Crow tephra at CRH 11 is at least MIS 7 in age (28). To confirm this, we collected four samples for additional OSL dating of this underlying unit.
OSL sampling focused on a 1-m-thick sequence of alternating point bar and overbank flow deposits interbedded with well-preserved macro-organic horizons. Samples CRH 11-1 and CRH 11-3 were both taken from 10- to 20-cm-thick, finely bedded sand lenses found immediately beneath the bone-rich gravels. CRH 11-2 was taken 60 cm below CRH 11-1, from a 20-cm sand lens underlying a dense organic mat of woody macros and leaves. CRH 11-4 was collected 40 cm above CRH 11-1 from within a 10- to 15-cm-thick, laterally continuous, gray clay-silt unit capping the main bone-bearing gravels. OSL samples were collected from cleaned exposures using opaque PVC tubes and wrapped in light-proof bags for transportation and storage. In the laboratory, quartz grains of 125–180 μm and 180–212 μm diameter were extracted from the unilluminated centers of the PVC tube samples under safe (dim red) light conditions and prepared for burial dose estimation using standard procedures (39), including etching by 48% hydrofluoric acid for 40 min to remove the α-irradiated external layers.
Single-grain equivalent dose (De) values were determined using the experimental apparatus described by Arnold et al. (40), the single-aliquot regenerative-dose (SAR) procedure (70) shown in Table S1, and the SAR quality assurance criteria outlined in Arnold et al. (41). A preheat of 160 °C for 10 s was used in the SAR procedure before measuring the natural, regenerative and test dose OSL signals. These preheating conditions yielded the most accurate measured-to-recovered dose ratios (1.01 ± 0.03) for ∼100 Gy dose-recovery experiments performed on multigrain (∼400-grain) aliquots of sample CRH 11-1 (Fig. S2A). Replicate single-grain dose-recovery tests performed on CRH 11-1 using a higher administered dose (∼315 Gy) yielded measured-to-recovered dose ratios consistent with unity (1.01 ± 0.03) (Fig. S2B) and relatively low levels of overdispersion (11 ± 6%), confirming the suitability of the chosen preheat conditions for high-dose single-grain De estimation.
Representative OSL decay curves and dose–response curves of grains that passed the SAR quality assurance criteria are shown in Fig. S3. Grains accepted for De analysis typically displayed sensitivity-corrected natural test-dose (Tn) signals of ∼200–300 counts/0.17 s (Fig. 3), although a small population of grains displayed Tn signal intensities of >500 counts/0.17 s (Fig. S3A). The majority of accepted grains display rapidly decaying OSL curves (reaching background levels within 0.4 s), which are consistent with quartz signals dominated by the most readily bleachable (so-called “fast”) OSL component. The single-grain dose–response curves are generally well represented by either a single saturating exponential function or a saturating exponential plus linear function, as has been widely reported for quartz grains with fast-dominated OSL signals (e.g., refs. 7173). A significant proportion of grains resemble the class I “supergrains” of Yoshida et al. (73) and display very high characteristic saturation dose (D0) values of 150–400 Gy (Fig. S3) when fitted with a saturating exponential dose–response function. These favorable dose–saturation properties have enabled us to obtain higher than average De values for these samples. Grains displaying continued linear growth at high doses also provided the opportunity to calculate finite De estimates with reasonable fitting uncertainties beyond the 2D0 practical limit of precise De interpolation suggested by Wintle and Murray (74).
Environmental dose rates were calculated using a combination of in situ FGS, HRGS, and ICP-MS, as detailed in Table S2. HRGS was performed on one of the samples (CRH 11-1) to assess the present-day state of secular (dis)equilibrium in the 238U and 232Th decay series of these sediments. The ratios of 226Ra/238U, 210Pb/226Ra, and 228Th/228Ra were 1.06 ± 0.05, 1.03 ± 0.06, and 0.98 ± 0.02, respectively, for this sample. These parent–daughter ratios are all within unity at 1σ or 2σ, indicating that a condition of secular equilibrium currently exists in the 238U and 232Th decay series.
The four samples share broadly similar single-grain De distribution characteristics (Fig. S4). The De datasets cover wide ranges of dose values (relative De range = 1.6–1.7), they are each consistent with a single dose population centered on the weighted mean De value (as indicated by the large proportions of grains lying within the gray bands in Fig. S4), and they are not significantly skewed at the 95% confidence interval when assessed using a log weighted skewness test (42, 75). These De distributions also display relatively low overdispersion values of 11–20% (Table S2), similar to the values of ≤20% reported for ideal, well-bleached single-grain samples that have not been affected by postdepositional mixing or β-dose spatial heterogeneity (e.g., ref. 42). These De distribution properties suggest that dose dispersion originating from extrinsic, field-related sources (e.g., partial bleaching, sediment mixing, β-dose heterogeneity) are relatively unimportant in relation to the size of the De measurement uncertainties. This interpretation is consistent with the sedimentological properties of these deposits (i.e., they comprise well-sorted sands and silts with clearly preserved primary sedimentary structures and boundaries). The broad range of De values obtained for the CRH 11 samples seems to be attributable to intrinsic, rather than external, sources of dose dispersion (i.e., scatter originating primarily from the experimental procedures themselves, particularly grain-to-grain variations in luminescence responses to the SAR conditions or the use of nonidentical field and laboratory bleaching/heating/irradiation conditions). This is borne out by the single-grain dose-recovery De dataset, which displays the same relative spread of values (dose-recovery relative De range = 1.6) as the natural De distributions and a consistent overdispersion value at 1σ (11 ± 6%). These similarities suggest that the specific sources of intrinsic De scatter captured by the dose-recovery test likely account for the natural De distribution characteristics of the CRH 11 samples.
In light of the aforementioned considerations, we have estimated the final burial doses of these samples from their weighted mean De estimates, calculated using the central age model (CAM) of Galbraith et al. (43). This choice of age model is also well-supported by the statistical criterion outlined by Arnold et al. (76); in all cases, the CAM yields a maximum log likelihood score (Lmax) that is statistically indistinguishable from that obtained using alternative age models with additional parameters [e.g., the minimum age models of Galbraith et al. (43) or the finite mixture model of Galbraith and Green (77)].
The single-grain OSL ages for the CRH 11 samples are in correct stratigraphic order and are statistically indistinguishable at 1σ (Table S2). The weighted mean OSL age for the lowest bone-bearing units (208 ± 6 ka) indicates sediment deposition during MIS 7 (78), and is consistent with the presence of the 124 ± 10-ka Old Crow tephra (29) in overlying deposits located ∼18 m above the present-day river level at CRH 11. The reliability of the MIS 7 OSL age for the lowest bone-bearing units at CRH 11 is also supported by the presence of prominent woody macro horizons within the dated fluvial sequence, which strongly indicate that interglacial rather than glacial conditions prevailed at the time of deposition.
Vertebrate fauna.
We visited CRH 11 over a 3-d period (June 27–30, 2008) to sample the lowest bone-bearing unit, consisting of a series of thin gravel beds separated by gray clays about 3.5 m above the then-current water level. We exposed this unit along a 10-m front at the base of the bluff, focusing on visible bone. Most were retrieved at the interface between one of the gravel layers, 1- to 2-cm thick, and an overlying gray clay. This and other evidence suggests that the depositional setting was an aggrading point bar, along which bones and teeth carried along by currents were dropped, likely near the proximal end of the point bar. Immediately below the bone-bearing stratum was an extensive bed of organic matter, featuring well-preserved leaves and twigs in dense mats. Small concentrations of organic matter also occurred as lenses within the gravel/clay sequence.
In this section we briefly compare the material we recovered with the comprehensive faunal lists for the Old Crow Basin (OCB) complied by Jopling et al. (27) in their table 2. Additional relevant notes on OCB fossil taxa may be found in Morlan (79) and references in Harington (80).
A total of 294 pieces recovered in the 2008 season were cataloged and incorporated into the collections of the Office of the Yukon Paleontologist. Virtually all of the pieces are fragments not readily assignable to species-level taxa. Isolated teeth were recovered, but these were rarely complete and often in very poor condition. This is consistent with the interpretation that the fossils were transported and rolled for a significant distance before final deposition on the point bar.
The list below consists of taxa positively recognized in this collection, together with Yukon Paleontology Program catalog numbers of representative elements (list is not exhaustive):
  • Mammuthus primigenius (319.0003, 0058, 0264: molar fragments). Although no complete teeth were found, the spacing between lamellae on some of the large molar fragments recovered suggests that M. primigenius is the mammoth species represented. Many of the numerous very large long bone fragments recovered are doubtlessly referable to this species as well.
  • Equus sp. (319.0015, 0068, 0277, 0280: molars and incisors). The taxonomy of the 40+ nominal species of Equus from Pleistocene North America is an unintelligible morass, all of the more so because ancient DNA evidence suggests that there may have been only one or two highly morphologically variable species in North America in the Late Pleistocene (81). In any case, no material diagnostic of any particular nominal species was recovered. Almost all horse teeth are split or otherwise damaged.
  • “Artiodactyl” (319.0257: astragalus). An artiodactyl astragalus too large to represent Rangifer may instead represent Soergelia or other ovibovin, although this cannot be affirmed in the absence of comparative material. This element has been sampled for ancient DNA purposes.
  • Rangifer tarandus (319.0282: premolar). Caribou is represented by isolated teeth and long bones.
  • Castoroides ohioensis (319.0189: incisor). Fragments of the distinctive, channeled incisors of this large beaver are well represented in the collection.
  • Castor canadensis (319.272: jaw). The much smaller living representative of Castor is definitively identified on the basis of a jaw fragment with distinctive p4 morphology.
  • “Rodent” (319.275; 319.271: incisor, tibia). A sliver of incisor enamel and a distal tibia, neither of which belong to either species of Castor, may represent at least one other big rodent (possibly Ochotona ?whartoni).
  • Gulo gulo (319.276: jaw). A left horizontal ramus with canine (broken) and P3 in situ. Measurements are close to those of extant wolverine and the fossil is therefore assigned to that species. This is the only definite carnivore element identified in the collection.
  • Megalonyx jeffersonii (319.013; jaw fragment). Although only the symphyseal area and an alveolus for the caninform is preserved on this fossil, it is adequate for species determination. (There is little else this specimen could represent except another megalonychid sloth, and none is known from arctic North America except M. jeffersonii.) This element has been sampled for ancient DNA purposes.
Jopling et al. (27) concentrated on developing the biostratigraphy, and possible archeological significance, of four OCB localities, designated as OCR (= CRH of this paper) 11, 12, 15, and 300. In total, the authors listed 35 mammal taxa at reasonably low taxonomic levels (species or genus), in addition to 5 taxa of birds and 6 of fishes, organized into 7 “faunal assemblages” (FAs). Two thirds of the mammal species listed are rodents. The remainder include a small number of carnivores, all still extant, and several large herbivores, of which Mammuthus, Equus verae, and Soergelia sp. are extinct.
The FAs were grouped and tentatively dated by Jopling et al. (27) in their table 2 as “Latest Illinoian,” “?Sangamonian,” or “Early Wisconsian,” depending on stratigraphic position and various age indicators. Although the total age range is reasonably large, the FAs are markedly oversplit, being largely defined by presence/absence of one or two supposed indicative species. Nevertheless, for convenience in reference the original FA designations are retained here.
The FAs of interest here, 1–4 in the Jopling et al. (27) scheme, are the oldest and occur in the stratigraphic layer denoted as Unit 1b at each of the four localities they studied. From section diagrams and descriptions provided in their paper, it is definite that our 2008 excavation took place in the same unit. FA 1 has the fewest listings, consisting of only three taxa at the species level, two of which are mammal (giant pika, Ochotona cf whartoni, and giant beaver, C. ohioensis). The lists for FAs 3 and 4 are also species-poor, but include Mammuthus, Equus, and Soergelia. FA 2 presents a much longer list, covering 26 positively or tentatively identified mammal taxa, including mammoth and horse but not steppe muskox. Bison was not identified in any of the early FAs, consistent with our finding as well. However, Jopling et al. did not have independent chronology for the site when they worked there because the Old Crow tephra had not been noted by those authors. Their first record of bison is in FA 6 is an important observation (noted as possibly Sangamonian by the authors). Those fossils would have been above the Old Crow tephra in our stratigraphy and thus of MIS 5e (Sanagamonian age) or later in age. We recovered no in situ fossils at these levels in our excavations.
M. jeffersonii, long known from fossils collected from various positions along the Old Crow River (82), was unaccountably omitted by Jopling et al. (table 2 in ref. 27) in their taxon lists. We confirm the occurrence of this taxon in Unit 1b at CRH 11 on the basis of an edentulous partial mandible. McDonald et al. (82) and Hoganson and McDonald (83) regarded most, and perhaps all, occurrences of Jefferson’s ground sloth in arctic North America as Sangamonian, but with the new dating of CRH 11 it is clear that this xenarthran had previously occupied this area during MIS 7.
Soergelia, the steppe muskox, is widely accepted as a Middle Pleistocene taxon. Although we did not find any material certainly referable to this taxon, Harington (84) was aware of its possible biostratigraphic significance as a potential Middle Pleistocene marker (“Kansan” in his terminology).

Ch’ijee’s Bluff.

Stratigraphy and chronology.
Ch’ijee’s Bluff (67° 29′ N, 139° 56′ W) is a prominent 4-km-long river exposure on the south bank of the Porcupine River, ∼15 km downstream of the village of Old Crow in northern Yukon. The bluff is within the continuous permafrost zone (>90% perennially frozen ground by area), although permafrost at the site thawed during the last glaciation when the region was inundated by large glacial lakes or tectonically controlled pluvial lakes (68, 69). There is also stratigraphic evidence for thaw of near-surface permafrost during MIS 5e, the last interglaciation (30).
The bluff exposes ∼50 m of unconsolidated Late Cenozoic sediment (Fig. 3) (31, 85). The lowermost third of the exposure comprises alluvial sand, silt, and gravel, with macrofossil and pollen biostratography suggesting a Late Pliocene age. Massive and crudely laminated silt and clay, likely associated with Early or Middle Pleistocene lacustrine sedimentation, make up the middle third of the bluff.
The upper third of exposed sediment at Ch’ijee’s Bluff hosts the unit of interest for this study. This unit, termed Unit 4 by earlier workers (31, 85), is composed of massive and crudely laminated silt with locally organic-rich interbeds and laminae. The Old Crow tephra (29) is present in this unit across much of Ch’ijee’s Bluff (Fig. S3A). Pollen, insects, and plant macrofossils from sediments directly associated with Old Crow tephra indicate shrub-tundra habitat and cooler than present climate (31). The tephra is consistently stratigraphically below a prominent bed of dark-brown, macrofossil-rich organic silt that commonly includes large spruce (Picea) logs. Pollen, insects, and plant macrofossils from this organic bed are indicative of closed boreal forest warmer than the present climate (31); spruce pollen percentage is higher than present-day surface samples, and the assemblage includes several extralimital taxa that are currently only found well to the south (e.g., the beetles Dyschirius laevifasciatus, Bradycellus lecontei, and Micropeplus sculptus; abundant bark beetles from the family Scolytidae; trace pollen from Typha; and macrofossils of the plants Carex sychnocephala, Chenopodium gigantospermum, and Sium suave). Logs and organic detritus from the organic bed all have nonfinite 14C ages (30, 31), indicating that they are older than the ∼50-ka limit of 14C dating and thus unlikely to be “young” organics placed into stratigraphic association with older sediments because of thaw slumping (86). Although Old Crow tephra is usually present 0.5–2 m below the organic bed, at some measured sections along Ch’ijee’s Bluff the organic bed truncates Old Crow tephra or grades into ice wedge casts that truncate the tephra, suggesting shallow ground thaw and associated sediment reworking (30). Collectively, at Ch’ijee’s Bluff and at other sites in eastern Beringia (38), the stratigraphic and paleoecological evidence indicate that Old Crow tephra was deposited late during the MIS 6 glacial, with the overlying organic bed representing warmer than present conditions during the last interglaciation (MIS 5e).
Bison at Ch’ijee’s Bluff.
In July 2006, we discovered a partial bison metacarpal at Ch’ijee’s Bluff (Fig. 3C) (YG 264.1) near the downstream limit of exposed Old Crow tephra at the bluff. The tephra was 40-m above river level during a period of relatively low water levels. At this measured section, the laterally continuous Old Crow tephra is 1- to 10-cm-thick and undulates with 20–30 cm of relief. The tephra is hosted in massive silt with local organic interbeds. A 5- to 15-cm-thick organic-rich bed with abundant logs and twigs sharply overlies the silts that host Old Crow tephra; this sharp contact is 1.5–2.5 m above the tephra. The bison metacarpal was found in situ in the silts, ∼125 cm above Old Crow tephra but several centimeters below the sharp contact with the woody organic-rich bed. Given the sharp contact between the tephra- and bone-bearing silt and the organic-rich bed, as well as the common occurrence of thaw-related unconformities associated with this sedimentary contact, the most parsimonious interpretation of the stratigraphic setting (Fig. 3B) is that the bison fossil dates to latest MIS 6, immediately preceding the MIS 5e interglaciation.

Reconstructing the Mitogenomic Evolutionary History of North American Bison.

To infer the evolutionary relationships among North American and Siberian bison, we performed phylogenetic analyses on a dataset of 46 bison and 4 yak mitochondrial genomes. This dataset comprised 10 present-day bison, one previously published Bison from Siberia, four previously published yaks, and 35 mitochondrial genomes that were isolated and assembled from Late Pleistocene and Holocene bison fossils as part of this study. Of these 35, 21 were not associated with either a radiocarbon date or stratigraphic information. We subsampled these and submitted them either to the KECK Accelerator Mass Spectrometry laboratory at the University of California (UC) Irvine or the Center for AMS (CAMS) at Lawrence Livermore National laboratory for ultrafiltered collagen extraction and radiocarbon dating. Sample details, ages, and GenBank accession numbers are provided in Dataset S1.

Mitochondrial genome reconstruction.

We performed all molecular biology work before the library indexing PCR step in the dedicated ancient DNA clean room at the UC Santa Cruz Paleogenomics laboratory, following strict protocols for ancient DNA research (87), for all bison other than the giant long-horned bison from a site near Snowmass, Colorado.
We extracted DNA from 100 to 250 mg of powdered bone using either of two silica-based methods suitable for the recovery of highly fragmented DNA (44, 45). DNA was extracted from 12 of the bison fossils during previous studies (Dataset S1 and references therein). We then prepared and amplified DNA libraries using the Meyer and Kircher (46) protocol, as modified by Heintzman et al. (88), from all 34 of these DNA extracts. We enriched these libraries for molecules that were similar to the plains bison (Bison bison) reference mitochondrial genome (GenBank RefSeq: NC_012346.1) via hybridization capture using custom, in-solution, biotinylated, 120mer or 80mer RNA MYbaits probes (MYcroarray). Enrichment experiments followed the manufacturer’s instructions (v1.3.8 or v2.3.1), with minimal modifications: the RNA baits were diluted 1:10 with 100 ng/μL yeast tRNA (Invitrogen) and we used the KAPA HiFi HotStart ReadyMix PCR kit (Kapa Biosystems) with a 30-cycle reaction, instead of the recommended 8–14 cycles, for the postenrichment amplification.
We sequenced the enriched libraries from both ends using either the Illumina MiSeq (UC Santa Cruz) or HiSeq-2000 (Vincent J. Coates Genomics Sequencing Laboratory at UC Berkeley) systems, for 2 × 75 or 2 × 150 cycles, respectively. We removed sequencing adapters and merged paired-end reads in SeqPrep (https://github.com/jstjohn/SeqPrep), with a minimum overlap of 10 bp required for reads to merge and discarded reads with a length shorter than 30 bp. We then mapped the merged and remaining unmerged reads to the B. bison reference mitochondrial genome using an iterative short-read assembler, MIA (50). To reduce the impact of DNA damage on the mitochondrial genome assembly, we only called consensus at each position where a minimum of three unique molecules mapped (3× coverage). If only three molecules mapped, all three had to be in agreement. Where more than three unique molecules mapped, we called a consensus at sites with at least 67% of called bases in agreement (51). Because mitochondria are haploid, this consensus-calling procedure, which takes into consideration the base quality score and the likelihood that a variant is the result of ancient DNA damage based on positional information in the molecule, provides an extremely robust means to determine consensus. Positions that did not meet these criteria were coded as missing data. This pipeline provided complete mitochondrial genome sequences for each library, with average coverage ranging from 6.6× to 898.4× (Dataset S1), with ∼75% of the samples having at least 20× average coverage.

The giant long-horned bison (Bison latifrons) from Snowmass, Colorado.

We obtained from the Denver Museum of Nature and Science a sample of the B. latifrons that was recovered from a recent excavation of the Ziegler Reservoir near Snowmass, Colorado (24). Snowmass is an ancient lake that formed at the top of a mountain ridge as glaciers retreated during the last part of MIS 6 (89). The chronology of the site has been established using a combination of radiocarbon and in situ cosmogenic 10Be and 26Al dates, uranium-series, and OSL dating (89, 90). We obtained a sample of a humerus of B. latifrons (DMNH EPV.67609), which was found in biostratigraphic zone 5d (BZ5d), and is associated with an age of 129 ± 10 to 113 ± 8 kyBP.
Because of its age and location of preservation, the Snowmass bison was, unsurprisingly, much more poorly preserved than the other bison fossils. We therefore used a highly sensitive approach to recover small fragments of DNA. In a dedicated ancient DNA clean room at the Max Planck Institute for Evolutionary Anthropology in Leipzig, Germany, we removed 50 mg of material using a sterile dentistry drill and extracted DNA following Dabney et al. (44). We then converted a fraction of the extract into a DNA library by single-stranded library preparation (47). We amplified and labeled the library with two sample-specific indices (91) and enriched for mitochondrial DNA in two rounds of bead-based hybridization capture (92), using a probe set that encompasses complete reference mitochondrial genome sequences from 242 eutherian mammal species, including B. bison (48). The enriched library was sequenced from both ends in a pool with other libraries on Illumina’s MiSeq using a 76+7+76+7-cycles recipe for double-index sequencing (91). Reads perfectly matching the expected index combination were isolated and full-length molecule sequences were reconstructed by overlap-merging of paired-end reads (93). Merged sequences longer than 30 bp were aligned to the plains bison reference mitochondrial genome using BWA (49) with “ancient” parameters (94). Duplicates were removed by calling a consensus from sequences with identical alignment start and end coordinates using bam-rmdup (https://bitbucket.org/ustenzel/biohazard). In total, we obtained 2,927 unique sequences that mapped to the plains bison mitochondrial genome (Dataset S1). To obtain a draft mitochondrial consensus sequence of the specimen, we followed a procedure used recently for consensus calling from highly degraded ancient DNA (51). Briefly, a consensus base was called for positions covered by at least three sequences that were at least in 67% agreement. We further masked T nucleotides in the first five positions of each sequence to reduce the impact of miscoding DNA damage on consensus quality. With this procedure, we obtained consensus calls for 10,714 positions of the specimen’s mitochondrial genome.

Assessment of mitochondrial DNA damage.

We scrutinized the reads mapped to bison for patterns of DNA damage, which are expected in degraded specimens (95). In a subset of samples, we assessed DNA fragment length distributions and the rates of miscoding lesions at the 5′ and 3′ ends of the aligned reads. These samples included the Ch’ijee’s Bluff steppe and Snowmass bison, as well as four other steppe bison from either Siberia (AE006, AE010) or North America (MS220, PH027). For the five steppe bison, we mapped the merged and remaining unmerged reads from SeqPrep to each steppe bison’s consensus sequence in BWA. This approach was to ensure that detected mismatches were not because of evolutionary divergence from the reference sequence. We analyzed damage patterns for each alignment using mapDamage v2.0.5 (52). For the Snowmass bison, we assessed DNA damage, using the unique reads aligned to the plains bison reference mitochondrial genome, in mapDamage.
The Ch’ijee’s Bluff and Snowmass bison samples display the shortest modal lengths (54 and 35 bp, respectively) (Fig. S6 A–F), consistent with their great age. As expected from depurination-induced fragmentation, all six samples exhibit noticeable increases in purine frequencies at the base position immediately preceding or following the 5′ and 3′ ends of reads, respectively (95), with the Snowmass bison sample most prominently displaying these patterns. Finally, all samples exhibit excesses of cytosine to thymine misincorporations at the ends of reads, respectively, consistent with the presence of deamination-induced base damage (95, 96). The Ch’ijee’s Bluff steppe bison sequences exhibit misincorporation rates similar to the other, younger steppe bison bone samples (Fig. S6 C, I, L, O, and R), whereas the Snowmass bison sequences display very strong accumulations of cytosine to thymine substitutions at their ends (70–73%) (Fig. S6F). We should caution that the damage pattern variability between samples is likely to be, at least in part, related to the diverse DNA extraction and library preparation methods used in this study.

Phylogenetic analysis.

We constructed an alignment of 46 bison mitochondrial genomes by supplementing our dataset (Dataset S1) with previously published mitochondrial genomes from 10 present-day North American bison [GU946976-84 (33); EU177871 (34)] and a ∼10,800-y-old bison from Rauchua River, Chukotka (KR350472) (35). As outgroup, we additionally included four yak (Bos grunniens; KJ463418, AY684273, KJ704989, KM233416) mitochondrial genome sequences, whose mitochondrial lineage is closer to ancient and modern North American bison than that of the European bison (Bison bonasus) (97). After aligning these sequences, which resulted in a 16,322-bp alignment, we created a second dataset in which we reduced the alignment to include only those positions where the base could be called with confidence for the Snowmass bison. This resulted in an alignment of 10,714 bp.
For both the full and reduced alignments, we extracted seven partitions containing the control region, 12S rRNA, 16S rRNA, concatenated tRNAs, and separate alignments for first, second, and third positions within the coding regions of the mitochondrial genome. We then used PartitionFinder v1.1.1 (53) to identify the simplest partitioning strategy that best represented the data and the best evolutionary model for each partition. For both alignments, PartitionFinder suggested four partitions: firstPositions+tRNAs+12S+16S (TrN+G for both alignments); secondPositions (HKY for the full alignment, HKY+I for the reduced alignment); thirdPositions (TrN+G for the full alignment, HKY+G for the reduced alignment), and the control region (HKY+I+G for both alignments).
We then performed genealogical inference using BEAST v1.8.3 (54), assuming the partitioning and evolutionary models described above. Because it allows sufficient flexibility to explore the coalescent history of the sequences included in this dataset, we assumed the skygrid coalescent prior (55) and the strict molecular clock (with separate rates for each partition), which we calibrated using age of each bison from which data were available. Ages of Late Pleistocene and Holocene bison were the median calibrated radiocarbon age before present (Dataset S1). Ages of the Ch’ijee’s Bluff and Snowmass bison were sampled from distributions (56). The age of the Ch’ijee’s Bluff steppe bison was sampled from a normal prior with a mean of 125 kyBP and SD of 4.5 kyBP, to reflect that it was found in situ between the Old Crow tephra and the MIS 5e forest bed. The Snowmass bison was sampled from a normal prior with a mean of 124 kyBP and SD of 8.5 kyBP, to reflect that it was found in a layer associated with MIS 5d. We ran two Markov chain-Monte Carlo chains for 60 million iterations each, sampling trees and priors every 3,000 iterations. Log and tree files were visually inspected using Tracer v1.6 (tree.bio.ed.ac.uk/software/tracer). The first 10% of sampled states were discarded as burn-in, after which the two runs were combined. Trees were summarized and the maximum clade credibility tree identified using TreeAnnotator, which is distributed as part of the BEAST package. The maximum clade credibility tree was visualized and annotated in FigTree v1.4.2 (tree.bio.ed.ac.uk/software/figtree).
The two alignments produced nearly identical phylogenies, with higher support values at the nodes when the reduced alignment was used (Fig. S5). The timing of bison entry into North America (node 1 in Fig. 1E) is inferred to be 193–138 kyBP (reduced dataset) or 172–134 kyBP (full dataset). The mitochondrial phylogeny can be roughly divided into three well-supported clades with strong phylogeographic partitioning. Siberian bison (green lineages) cluster together in the most deeply diverging clade, within which we identify a second, strongly supported wave of bison dispersal into North America (node 2 in Fig. 1E). This second dispersal occurred during the period 45–21 kyBP (reduced dataset) or 36–21 kyBP (full dataset), broadly within MIS 3 or 2 (Fig. S5), during which the Bering Land Bridge was exposed (36). We note that although the second wave clade is poorly sampled in our dataset, and the timing may change slightly with additional samples, two waves of dispersal are clearly indicated by these data. Extant American bison mitochondrial diversity seems to be derived exclusively from the first wave of bison dispersal into North America (Fig. 1E and Fig. S5). We further note the single deeply diverging lineage from Rauchua River, Siberia, which was recently isolated from a well-preserved partial mummy dating to ∼10.8 kyBP (35). Interestingly, this is the youngest Asian steppe bison from which mitochondrial data are available, and suggests that only a small subset of Asian diversity ever crossed the Bering Isthmus into North America. We do not find evidence for the migration of North American bison back into Asia, which tentatively suggests unidirectionality in Pleistocene bison migration across the Bering Land Bridge.

Materials and Methods

This section provides an overview of the methods of this study; full details can be found in SI Text.

Geochronology.

Chronology at Ch’ijee’s Bluff and CRH 11 relies on identification of tephra in sediment exposures and optically stimulated luminescence (OSL) dating. Old Crow tephra was identified based on stratigraphic position, glass shard morphology, and grain-discrete glass major element geochemistry. Glass geochemical analyses were by wavelength dispersive spectrometry on a JEOL 8900 electron microprobe at University of Alberta following Reyes et al. (38), with correlations confirmed by concurrent analyses of an Old Crow tephra reference sample (Fig. S1).
We obtained four samples for single-grain OSL dating from immediately above and below the lower fossil-bearing horizon at CRH 11. Samples were processed under safe (dim red) light conditions using standard procedures (39) to isolate refined quartz fractions. We performed equivalent dose (De) measurements on 1,800–2,400 individual quartz grains per sample using the experimental apparatus described by Arnold et al. (40) and the single-aliquot regenerative-dose procedure shown in Table S1. We considered 3–6% of the measured grain populations suitable for De determination after applying the SAR quality assurance criteria (41). De estimation over high dose ranges (>300 Gy) was well-supported by the single-grain dose saturation characteristics and dose-recovery test results (Figs. S2 and S3). The natural De datasets exhibited low overdispersion values (11–14%) and are dominated by experimental rather than field-related De scatter (42) (Fig. S4). We therefore calculated the final burial doses using the central age model (43). Dose rates were calculated using a combination of field γ-ray spectrometry (FGS), high-resolution γ-ray spectrometry (HRGS), and inductively coupled plasma-mass spectrometry (ICP-MS) (Table S2). To calculate the final OSL ages, we assumed that the measured radionuclide activities and present-day field water/organic contents prevailed throughout the burial period of these perennially frozen deposits. An uncertainty of 10% was assigned to long-term water and organic content estimates to accommodate minor variations during the burial periods.
Table S1.
SAR procedure used in this study
Table S2.
Dose rate data, equivalent doses (De), and OSL ages for the CRH 11 samples

DNA Extraction, Sequencing, and Mitochondrial Genome Assembly.

We assembled mitochondrial genomes for 35 ancient bison, including the Ch’ijee’s Bluff steppe bison and a giant long-horned bison from a site near Snowmass, Colorado (24). Of these, 21 were not associated with any stratigraphic or age information, and were sent to accelerator mass spectrometry (AMS) radiocarbon dating facilities for dating using ultrafiltered collagen (Dataset S1). We extracted DNA from 23 ancient bison using silica-based methods optimized for recovery of ancient DNA (44, 45), and included in our dataset 12 previously extracted bison (Dataset S1, and references therein). We converted extracted DNA to either double-stranded (46) or single-stranded (47) Illumina-compatible libraries. Mitochondrial DNA molecules were enriched using biotinylated RNA baits based on either the bison mitochondrial genome (GenBank: NC_012346) or a 242 mammal mitochondrial genome reference panel (48). We sequenced enriched libraries on the Illumina MiSeq or HiSeq platforms using paired-end chemistry. Sequencing read pairs were merged and adapter trimmed in SeqPrep. Merged and remaining unmerged reads were mapped to the bison mitochondrial genome using either Burrows–Wheeler Aligner (BWA) (49) or the iterative short-read assembler, MIA (50). We collapsed PCR duplicates using either bam-rmdup or MIA. For consensus sequence calling, we required each position to have a minimum of 3× coverage and a base agreement greater than 67% (51). To evaluate DNA preservation in these oldest bison, we used mapDamage (52) to assess patterns of DNA fragmentation and cytosine deamination (Fig. S6). The resulting ancient mitochondrial genomes ranged in coverage from 6.6× to 898.4× (Dataset S1). The Ch’ijee’s Bluff steppe bison had an average fragment length of 54 bp and was sequenced to 159× coverage, with all bases called following consensus calling, as above. The Snowmass bison, which was much more poorly preserved (Fig. S6 D–F), had an average fragment length of 35 bp, with 5,596 missing bases following consensus calling.
Fig. S6.
DNA damage patterns for the Ch’ijee’s Bluff steppe bison (MS226/DC009) (A–C), B. latifrons from Snowmass (B5493) (D–F), two Siberian steppe bison (AE006, G–I) (AE010, J–L), and two steppe bison from North America (PH027, M–O) (MS220, P–R). Damage patterns include DNA fragment length distributions (A, D, G, J, M, P), fragmentation patterns (B, E, H, K, N, Q), and misincorporation rates (C, F, I, L, O, R), the latter of which are cytosine to thymine (red lines) and guanine to adenine (blue lines) misincorporations. PH027 was previously reported in Heintzman et al. (98).

Phylogenetic Analysis.

We aligned mitochondrial genomes from the 35 ancient bison described above with previously published mitochondrial genomes from 10 present-day American bison (33, 34), an ancient Siberian bison (35), and 4 yaks (Bos grunniens). We then created two datasets for phylogenetic analysis: one comprising the complete mitochondrial genome (full dataset), and another limited to only those sites in the mitochondrial genome in which we were able to call a consensus base for the Snowmass bison (reduced dataset). We partitioned both alignments and selected appropriate models of molecular evolution using PartitionFinder (53), and inferred the evolutionary relationships among the sampled mitochondrial lineages using BEAST v1.8.3 (54). Following model testing (SI Text), our final analyses assumed a strict molecular clock and the skygrid coalescent prior (55). We calibrated the molecular clock using median calibrated radiocarbon ages for each sampled mitochondrial genome, and sampled the ages of Ch’ijee’s Bluff and Snowmass bison using a mean and SD of 125 ± 4.5 kyBP and 124 ± 8.5 kyBP, respectively (56). For each analysis, we ran two Markov chain-Monte Carlo chains for 60 million iterations each, sampling priors and trees every 3,000 iterations, and discarding the first 10% as burn-in, and combining the remainder. We visually inspected log files for run convergence using Tracer and summarized the sampled trees using TreeAnnotator. Phylogenies presented in the text and SI Text are maximum clade credibility trees (Fig. 1E and Fig. S5).

Data Accessibility.

Dataset S1 includes repository and radiocarbon accession details for all fossil specimens analyzed. The Ch’ijee’s Bluff steppe bison specimen, YG 264.1, is archived in the fossils collections of the Vuntut Gwitchin First Nation Government in Old Crow, Yukon. The giant long-horned bison specimen from the site near Snowmass, DMNH EPV.67609, is archived in the vertebrate paleontology collections of the Denver Museum of Natural History, Colorado. Mitochondrial genome sequences have been deposited in GenBank, with accession numbers KX269109, and KX269112-KX269145. Input BEAST files are available as Dataset S2.

Acknowledgments

We thank Andrew Fields, Brandon Letts, Dan Chang, Joshua Kapp, and Amanda Hart for technical assistance; Ben J. Novak for the bison skull paintings in Fig. 2; Ian Miller, Joe Sertich, and the Denver Museum of Nature and Science for access to the Snowmass sample; and Harold Frost, Jr. and the late Stephen Charlie for assistance in the field. Access and support to work at the Old Crow sites was facilitated by the Heritage Department of the Vuntut Gwitchin Government. This work used the Vincent J. Coates Genomics Sequencing Laboratory at University of California, Berkeley, supported by NIH S10 Instrumentation Grants S10RR029668 and S10RR027303. Support was provided from Yukon Heritage Branch; Natural Resources Canada Polar Continental Shelf Program (D.F.), Canada Research Chairs program, and Natural Sciences and Engineering Research Council (D.F.); National Science Foundation ARC-09090456 and ARC-1417046, and Gordon and Betty Moore Foundation GBMF3804 (to B.S.); and Australian Research Council DP0878604 and FT130100195 (to L.J.A.).

Footnotes

  • Author contributions: D.F. and B.S. designed research; D.F., M.S., P.D.H., A.V.R., G.D.Z., A.E.R.S., M.M., E.H., B.J.L.J., L.J.A., R.D.E.M., and B.S. performed research; G.D.Z. and L.J.A. contributed new reagents/analytic tools; P.D.H., A.V.R., A.E.R.S., and L.J.A. analyzed data; and D.F., P.D.H., A.V.R., G.D.Z., M.M., L.J.A., R.D.E.M., and B.S. wrote the paper.
  • The authors declare no conflict of interest.
  • This article is a PNAS Direct Submission.
  • Data deposition: The sequences reported in this paper have been deposited in the GenBank database (accession nos. KX269109, and KX269112KX269145).
  • This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1620754114/-/DCSupplemental.

References

View Abstract

Nenhum comentário:

Postar um comentário

Observação: somente um membro deste blog pode postar um comentário.