During his expeditions to the North Pole Robert Peary harvested at least 233 Peary caribou (Rangifer tarandus pearyi) from the northeast (NE) portion of Ellesmere Island, Nunavut, Canada. To assess the impact of this intensive local harvest we compared historic and contemporary samples. We explored 2 models of repopulation for this region. In the recovery model animals not hunted by Peary were the source of the current population. In the recolonization model animals from adjacent regions contributed to the current population. A genetic assessment of mitochondrial DNA (mtDNA) diversity comparing historic specimens (n = 12) to contemporary samples from the same region (n = 22) indicates the loss of 1 major mtDNA haplogroup. In this region historic and contemporary samples were significantly differentiated from one another. Combined mtDNA (n = 121) and 9 nuclear DNA loci (n = 151) indicate that the contemporary NE animals are most similar to those in adjacent regions to the south. We interpret this as supporting our recolonization model and suggest that, following Peary, animals dispersed north from central Ellesmere. Animals from the region to the west (NW) of NE showed the greatest differentiation from all other regions, suggesting that although proximate in distance, movement between the NW and NE might always have been limited. The NW region was unique in that it contained 1 predominant mtDNA haplotype at high frequency (0.88), the lowest observed heterozygosity (HO = 0.50), and the highest mean relatedness (R = 0.063). Our combined results highlight the sensitivity of populations at the range limit to stochastic events and the potential limits to recovery following significant perturbations.
Genetic implications of a severe population reduction have been well characterized, both theoretically (Lacy 1997; Lande and Shannon 1996) and empirically (Lacy 1987; Pastor et al. 2004). Evidence suggests that loss of genetic variability in a population can lead to lower fitness, especially in the face of environmental stress (Coltman et al. 1999; Kristensen and Sørensen 2005; Woodworth et al. 2002). Confounding these expectations are examples of species with currently large populations that have recovered from a severe bottleneck but show little variation. Examples include northern elephant seals (Mirounga angustirostris—Hoelzel 1999) and moose (Alces alces—Broders et al. 1999). Examples also exist of species that have persisted for long periods with low genetic diversity such as Kodiak Island brown bears (Ursus arctos—Paetkau et al. 1998). Other populations have experienced a demographic bottleneck with no detectable reduction in genetic diversity; examples include banner-tailed kangaroo rats (Dipodomys spectabilis—Busch et al. 2007) and wild rabbits (Oryctolagus cuniculus—Queney et al. 2000).
Large mammals may be at greater risk than many other groups to experience bottleneck events because they tend to have lower population sizes and lower rates of increase, in addition to being resources or competitors of many human populations. Long generation times characteristic of large mammals can extend the time needed for a population to recover, but they also act to maintain genetic diversity through a bottleneck by reducing the number of generations that experience the bottleneck. In addition, traits such as high dispersal ability and overlapping generations can act to increase population-wide effective size (Ne) and buffer against recurrent population bottlenecks (Busch et al. 2007; Kuo and Janzen 2003).
To further examine questions of bottlenecks and potential recovery of a large mammal species we used molecular genetic techniques to examine the recent (100 years to present) demographic history of a large ungulate, Peary caribou (Rangifer tarandus pearyi), at its northern range limit. A minimum of 233 Peary caribou from northern Ellesmere Island, Canada, were harvested by Robert E. Peary to support his 2nd (1905–1906) and 3rd (1908–1909) North Pole expeditions, with 149 caribou being harvested during the final expedition (Dick 2001; Manseau et al. 2004). Samples are available from this harvest through the American Museum of Natural History, New York.
Peary caribou occur primarily on islands throughout the Canadian Arctic Archipelago north of the 74th parallel. They have experienced declines in population size over much, if not all, of their range during the last 4 decades (Committee on the Status of Endangered Wildlife in Canada [COSEWIC] 2004) from an estimated high of 28,288 (SE ± 2,205) in 1961 (Miller et al. 2004) to an estimated low of 2,500 in the late 1990s (COSEWIC 2004). Significant declines were reported over the winter of 1973–1974 (Miller and Russell 1975; Parker et al. 1975) and in the late 1990s (Miller and Gunn 2003; Miller et al. 2007), and were associated with winter rain events and deep snow that prevented caribou from accessing forage (Miller and Gunn 2003; Parker et al. 1975). These declines were best documented on western islands (i.e., Banks Island and Bathurst Island) and assumed to have occurred on all islands (COSEWIC 2004), although some experts have commented that the Queen Elizabeth Islands in the eastern Arctic and especially Ellesmere Island likely have had relatively stable population numbers (Taylor 2005). Other factors leading to the declines, including overharvesting and competition with muskox, also have been suggested (COSEWIC 2004; Vincent and Gunn 1981) but likely play only a minor role in most areas (but see Miller et al. 2007). The overall population reduction led the Committee on the Status of Endangered Wildlife in Canada (COSEWIC) to list Peary caribou as Threatened in 1979 (Gunn et al. 1979) and subsequently as Endangered (COSEWIC 2004; Miller 1991). Peary caribou also are listed as Endangered by the World Conservation Union (International Union for the Conservation of Nature and Natural Resources [IUCN] 1996).
Ellesmere Island is the largest and northernmost island in the Canadian Arctic Archipelago. Ellesmere is the northern part of the range for Peary caribou and because of its remote location, limited information on population size and trends is available. The most recent survey conducted in 2005–2006 counted 587 animals (Campbell 2006). An estimated population size is not yet available, but this clearly points, as did previous survey efforts, to low densities (Gauthier 1996; Manseau et al. 2005; Riewe 1973) and clustered distribution throughout the island. Records of significant harvest activities are available for northeast (NE) Ellesmere Island and these harvests consist primarily of those of early European explorers (Manseau et al. 2004, 2005). Between 1898 and 1910 Robert Peary's expeditions harvested at least 233 caribou from northern Ellesmere Island between Lake Hazen and the northern coast (Dick 2001; Manseau et al. 2004; Fig. 1). Between 1910 and the present only 1 hunting record exists for northern Ellesmere Island; 3 caribou were killed in 1935 as part of the Oxford University Ellesmere Island Expedition (Manseau et al. 2005).
Figure 1
Map of northern Ellesmere Island with track line of Robert Peary's hunting expeditions and caribou kill sites (Manseau et al. 2004). Historic samples used in this study are from kill sites in the northeast portion of Ellesmere Island.
Northeast Ellesmere Island contains Quttinirpaaq National Park and the military enclave of Alert. Surveys by Parks Canada between 1988 and 2002, observations by park wardens, and radiocollaring work in the late 1990s within Quttinirpaaq National Park suggested a minimum estimate of 85 animals for the northern part of Ellesmere Island (Manseau et al. 2004). By contrasting the number of caribou that Peary harvested with recent estimates we infer that NE experienced a significant population reduction at the turn of the 20th century and that the animals currently found in this region are either descendents of caribou that survived the Peary harvests (the recovery model) or animals that have dispersed in from adjacent regions (the recolonization model). Under the recovery model we predict that genetic diversity will be low compared to other regions, and the signature of a population bottleneck will be observed. Additionally, low Ne during recovery could result in increased differentiation from adjacent areas due to genetic drift. Under the recolonization model we predict genetic diversity will be comparable to other regions, and specifically, the current sample will be genetically similar to a neighboring (source) population.
Samples from the NE region were compared with the northwest (NW), central (CE), and south-central (SC) regions. The NW region encompasses the Marvin Peninsula and is separated from the NE, CE, and SC regions by mountainous and glaciated terrain (Fig. 2). No record exists of exploration or of hunting in the NW region, and therefore these animals have probably existed unperturbed for generations. The CE and SC regions have generally higher productivity, based on AVHRR-NDVI satellite imagery (Manseau et al. 2004), which may contribute to the higher density of caribou. Although some Peary caribou have been harvested in these regions over the last 100 years, a potential population bottleneck would likely be caused by extreme weather events as reported for southern islands.
Figure 2
Locations of Peary caribou samples collected in 2004 and 2006. Four regions discussed in the text are northeast (NE, area north of the southern end of Lake Hazen, not including NW), northwest (NW, Marvin Peninsula near Ward Hunt Island), central (CE, area between Nansen Sound and Greely Fiord to the southern end of Lake Hazen), and south-central (SC, area south of Nansen Sound and Greely Fiord).
The opportunity to assess the impact of a highly localized but significant harvest on the NE population is enhanced further by the analysis of historic hides collected by Peary's expedition. The use of temporal replicates to assess genetic changes in populations has increased as methods for extracting DNA from historic and ancient sources have improved (Leonard et al. 2005; Martínková and Searle 2006). These methods typically use mitochondrial DNA (mtDNA) sequences to examine haplotype frequencies before and after a documented bottleneck event. We hypothesize the loss of haplotype diversity under the recolonization model or changes in haplotype frequencies under the recovery model. These relatively simple models could be confounded by the intervening 100 years during which any number of population crashes, expansions, extirpations, or recolonization events might have occurred. Fortunately, comparison of caribou from NE to those from adjacent regions to the west and to the south is possible.
Materials and Methods
Sample sources
We obtained skin samples (tanned hides) from the American Museum of Natural History, New York. These samples were originally collected during the Peary expeditions in 1905 (n = 42) and 1908 (n = 6). Location data for each specimen were determined from original labels or from journals kept by expedition members (Manseau et al. 2004; Fig. 1). These samples will be referred to as the historic samples.
Fecal pellets collected in 2004 and 2006 by Parks Canada and the Government of Nunavut provided the source DNA for the contemporary samples. When these samples were collected geospatial positions (latitude and longitude) were recorded (Fig. 2). Fecal pellets were sampled during population surveys conducted using Bell 206L helicopters (Bell Helicopter Textron Inc., Hurst, Texas). In August 2004 systematic transect lines were flown over northern Ellesmere Island (Manseau et al. 2005). In April–May 2005 and 2006 systematic transects were surveyed over most of the unglaciated land area of Ellesmere Island covering an estimated 40,000–50,000 km of transect lines (Campbell 2006).
Historic samples
From each historic sample we used a 2-cm2 piece of hide to obtain DNA. The hide was rehydrated in 1× phosphate-buffered saline (pH 7.4) then moved to lysis buffer (4 M urea, 0.2 M NaCl, 0.5% n-lauroyl sarcosine, 10 mM 1,2-cyclohexanediaminetetraacetic acid, and 0.1 M Tris HCl pH 8). Samples were digested using sequential additions of Proteinase K (20 mg/L; Roche, Indianapolis, Indiana) in the following order; 15 units (U) followed by incubation for 1 h at 65°C; 15 U followed by 24 h at 37°C; additional Proteinase K was added to samples that were not digesting completely.
During all steps we took precautions to prevent and identify contamination. An isolated fume hood was used for the initial sample preparation. This area was thoroughly cleaned with Decon 75 (Fisher Scientific, Pittsburgh, Pennsylvania) and bleach before and after each use. Samples then were transferred for extraction to a paleo-DNA laboratory at Trent University; this laboratory had no previous exposure to caribou DNA. We extracted samples in small batches (6–12), and extraction negatives were carried through for each batch.
We explored several extraction strategies; however, the most effective (resulting in the greatest number of successful amplifications) consisted of 2 sequential extractions. Initially, we isolated DNA from a large volume of digested hide (2 ml) using traditional phenol ∶ chloroform extraction. Amplification of this extract tended to fail due to the presence of polymerase chain reaction inhibitors; therefore, we employed a 2nd extraction of the DNA from the phenol ∶ chloroform extraction using a DNeasy tissue kit (Qiagen Inc., Valencia, California) to remove inhibitors (Calvignac et al. 2008).
Amplification reactions were set up in the non-DNA module of the paleo-DNA laboratory using reagents dedicated to this project. Sealed reactions then were transferred out of the paleo-DNA laboratory for polymerase chain reaction amplification, after which amplified product did not return to the paleo-DNA laboratory. We amplified all extraction negatives with each primer set to identify contamination, and only batches with clean extraction negatives were carried forward.
We designed 3 primer sets to amplify overlapping fragments of the mitochondrial control region; each amplicon was approximately 200 base pairs (bp) in length (Table 1). Amplification of each fragment was carried out in a total volume of 40 µl. Final concentrations of the following ingredients were used: 1× PCR buffer (Invitrogen, Carlsbad, California), 1.5 mM of MgCl2, 0.2 mM of deoxynucleoside triphosphates, 0.4% bovine serum albumin, 0.2 µM of each primer, and 1 U of Taq polymerase (Invitrogen). Four microliters of template DNA was used per reaction. Thermocycling conditions were as follows: 94°C for 5 min, 35 times (94°C for 60 s, 50°C for 60 s, 72°C for 60 s), and 72°C for 2 min. If amplification resulted in insufficient product for sequencing, a subsequent reaction was conducted using DNA isolated and gel-extracted from the 1st reaction. Concentrations and reaction conditions for the reamplification were modified to reduce bovine serum albumin to 0.1% and the number of cycles to 29. Amplified products were cleaned using Exosap IT (Invitrogen) and sequenced in both directions using a MegaBace 1000 capillary sequencer (GE Healthcare, Piscataway, New Jersey).
Table 1
Primers used to amplify control region in historic and recent samples of Peary caribou.
Contemporary samples
Fecal samples from the 2004 collection (n = 56) consisted of pellets that had been washed in 2 ml of phosphate-buffered saline, whereas the 2006 collection (n = 156) consisted of frozen fecal pellets. Extraction volumes differed to account for increased initial volume of the washed samples. We used the Qiagen DNAeasy tissue kit to recover DNA from epithelial cells that were suspended in the wash buffer and from swabs taken from the frozen pellets (Ball et al. 2007). We quantified our extracted double-stranded DNA using the BMG FluoStar Galaxy (BMG LabTech Inc., Cary, North Carolina) fluorescence absorbency PicoGreen assay (Molecular Probes Inc., Eugene, Oregon).
We amplified a portion of the mitochondrial control region in contemporary samples using primers L15693 and H00068 (Hundertmark et al. 2002). Reagent concentrations were as outlined above in a 25 µl total volume and using approximately 5 ng of DNA template. Thermocycling conditions were as above but modified to perform fewer (29) cycles. Amplification products were cleaned using Exosap IT (Invitrogen) and sequenced in both directions using a MegaBace 1000 capillary sequencer (GE Healthcare). Sequences were edited and aligned by eye in the alignment module of MEGA3.1 (Kumar et al. 2004).
We generated genotype profiles for all contemporary samples using 11 microsatellite loci: RT5, RT6, RT7, RT9, RT24, and RT30 (Wilson et al. 1997), and Map2C, BL42, BMS1788, BM888, and BM848 (Bishop et al. 1994). We employed a multiplexing and pooling strategy to obtain genotypes (see Appendix I for details, including primer concentrations and annealing temperatures). Pooled samples were profiled using a MegaBace 1000 capillary sequencer (GE Healthcare). Multiplex reactions contained concentrations of ingredients in a total volume of 10 µl: 1× PCR buffer (Invitrogen), 2 mM of MgCl2, 0.2 mM of deoxynucleoside triphosphates, 0.2% bovine serum albumin, and 1 U of Taq polymerase (Invitrogen). Approximately 5 ng of template DNA was used per reaction. Thermocycling conditions were as follows: 94°C for 5 min, 29 times (94°C for 30 s, Tanneal for 60 s, 72°C for 60 s), and 60°C for 45 min. Samples were sexed using amplification of the zinc finger gene with the primers SDP730 (5′-GGA AAT CAT TCA TGA ATA TCA C-3′—Ball et al. 2007) and LGL 335 (5′-AGA CCT GAT TCC AGA CAG TAC CA-3′—Shaw et al. 2003).
Levels of genetic diversity
Mitochondrial diversity was measured using haplotype frequencies (Nei 1987) and nucleotide diversity in the program Arlequin version 3.11 (Excoffier et al. 2005). The model of sequence evolution used to calculate nucleotide diversity was obtained using MODELTEST version 3.7 (Posada and Crandall 1998) based only on the larger contemporary control region fragment. The best model found using the Akaike information criterion (AIC—Akaike 1981) was the K81uf +I model (AIC = 1,286), which includes parameters for variable base frequencies, different transition and transversion rates, and the proportion of invariable sites. Because this option is not implemented in Arlequin 3.11, the Tamura–Nei model (Tamura and Nei 1993), for which the AIC score was essentially equivalent (AIC = 1,288), was used. This model also considers base frequencies and substitution rates to correct for multiple nucleotide substitutions at the same site.
Because individual caribou can be sampled multiple times when using noninvasive fecal sampling, individuals were identified using genetic profiles and the program GENECAP (Wilberg and Dreher 2004). Only 1 sample from each individual was used for subsequent analyses. Mean number of alleles, expected unbiased heterozygosity (Nei 1987), and observed heterozygosity were calculated with the Excel Microsatellite toolkit (Park 2001). Allelic richness, which standardizes the number of alleles to be independent of sample size, was calculated with FSTAT version 2.9.3.2 (Goudet 1995). GENEPOP version 4.0 (Rousset 2008) was used to assess linkage disequilibrium and deviations from Hardy–Weinberg equilibrium. Sequential Bonferroni corrections (Rice 1989; Rousset 2008) were conducted to account for multiple tests.
We assessed the similarity of individuals within the entire data set and within each region using 3 measures; FIS, relatedness, and probability of identity. FSTAT version 2.9.3.2 was used to calculate FIS, and GENALEX version 6.1 (Peakall and Smouse 2006) was used to calculate relatedness and probability of identity. Relatedness was calculated using the formula of Ritland (1996) because this measure is appropriate when few loci are used and when loci have <6 alleles (Ritland 2000). Permutation (999) and bootstrap (1,000) replications were conducted to identify the 95% confidence intervals around the null hypothesis of relatedness = 0 and the observed mean relatedness for each region. Probability of identity was calculated assuming a random association of loci and alleles (pI) and the more conservative probability of identity for siblings (pIsibs—Waits et al. 2001).
Population differentiation among regions was assessed using FST (Weir and Cockerham 1984) and the program FSTAT version 2.9.3.2. Significance was tested using 1,000 permutations and inferred after Bonferroni correction to account for multiple tests.
Historic NE to contemporary differences
Mitochondrial control region differentiation between temporal replicates (historic versus contemporary NE) and among contemporary regions was calculated using pairwise FST measures in the program Arlequin version 3.11. We interpreted significant differences between historic NE and contemporary samples as indicating divergence since 1908. Nonsignificant pairwise values among contemporary sampling areas were used to infer the current patterns of gene flow among caribou. Population differentiation was determined using pairwise FST-values calculated from mtDNA haplotype frequencies alone and from the distance matrix constructed using the Tamura–Nei algorithm (Tamura and Nei 1993). We assessed significance with 1,000 permutations of the data.
Population bottlenecks
Population bottlenecks can be detected in genetic data because at low population numbers the effect of genetic drift is more pronounced. We employed 2 methods to detect bottleneck events in the recent past of Peary caribou on Ellesmere Island. Evidence of a recent bottleneck was tested using the program BOTTLENECK version 1.2.02 (Cornuet and Luikart 1996; Piry et al. 1999). Generally, when a population has experienced a significant reduction in population size it will lose rare alleles. However, these alleles do not contribute significantly to measures of heterozygosity and a bottleneck will cause an excess of heterozygosity over what would be expected given the number of alleles observed (Cornuet and Luikart 1996). BOTTLENECK compares the expected heterozygosity, given the number of alleles observed, to the observed heterozygosity and tests this difference using a Wilcoxon signed-rank test (Piry et al. 1999). Identification of a bottleneck is sensitive to which underlying mutation model is assumed; 3 models were explored that represent liberal (infinite alleles model), intermediate (2-phase model), and conservative (stepwise mutation model) alternatives (Cornuet and Luikart 1996). We implemented the 2-phase model assuming 95% stepwise mutations as suggested by Piry et al. (1999). We performed these tests for each region and with all regions combined.
The M-ratio proposed by Garza and Williamson (2001) also was used to detect the signature of a population bottleneck. This statistic identifies the signature of a population bottleneck by examining the allele frequency distribution for larger-than-expected gaps in the allele range, given the number of alleles (Garza and Williamson 2001). If the population size remains low after a bottleneck, Garza and Williamson (2001) suggest that the signature of a bottleneck will be detectable longer using the M-ratio than using BOTTLENECK. Based on empirical evidence, M-ratio values less than 0.7 are generally indicative of populations that have experienced a recent reduction in size (Garza and Williamson 2001). The M-ratio can be tested for significance by comparing the observed value to a simulated population at equilibrium. This requires several parameters to be specified. The percentage of mutations greater than a single step was set at 10%, and the mean step size for those mutations was set at 3.5, as suggested by Garza and Williamson (2001). The program M_P_Val (Garza and Williamson 2001) was used to calculate the M-ratio.
Results
Data characteristics
Twelve (21%) historic NE samples (AMNH accession numbers in Appendix II) were amplified successfully and sequenced for a 214-nucleotide mtDNA fragment with no missing or ambiguous bases. This fragment was amplified with primer sets 1 and 2 (Table 1). Successful amplification of 121 (80%) contemporary samples resulted in a 427-nucleotide fragment. However, for consistency between the historic NE and contemporary data sets, all analyses were conducted using only the 214-bp fragment. This mtDNA data set consisted of 133 samples and contained 11 variable sites (Table 2).
Table 2
Definitions for the haplotypes identified in a 214-bp fragment of control region sequence. Number of individuals with the observed haplotype are reported for the historic NE sample, each contemporary region, and the combined contemporary sample (All). Regions of Ellesmere Island are abbreviated as follows; northeast (NE), northwest (NW), central (CE), and south-central (SC).
Two hundred twelve fecal samples were genotyped at 11 microsatellite loci and represented 151 individuals. The sex ratio was even in the total data set (proportion of females = 0.48), with a slight male bias in NE (proportion of females = 0.32; Table 3). Two loci were dropped from further analysis because they essentially were fixed (BM888; 80% allele 179) or were missing data in many samples (RT7; not profiled in 58%). Significant heterozygote deficits were identified in 4 loci after Bonferroni correction (RT6, RT9, RT24, and Map2C). When each region was considered independently only RT9 (NE and SC) and Map2C (SC) had significant heterozygote deficits. Further investigation of these deviations using the program MICRO-CHECKER (van Oosterhout et al. 2004) identified only RT24 and RT6 as having a significant heterozygote deficits. These deviations can be a result of biological characteristics that violate the assumptions of Hardy–Weinberg equilibrium, for example nonrandom mating. Alternatively, they could be the result of mutations in the priming region that cause the failure of alleles to amplify (null alleles) or problems associated with alleles failing to amplify due to technical issues (drop out). Null alleles are not implicated with these loci because Zittlau (2004) used RT24 in multiple herds of Peary caribou and did not observe heterozygote deficiency, although she did find evidence for heterozygote deficiency at RT9 in 1 population. Allelic drop out seems unlikely, because multiple extractions and amplifications of samples generally did not alter the genotype. Therefore, these 2 loci were kept in the analysis. Obviously, these deviations will affect statistics that assume Hardy–Weinberg equilibrium. In these cases calculations were done with and without RT9 and RT24 and reported when this action caused a change in our results.
Table 3
Sample size (n), proportion of females in sample (pF), mean number of alleles (Na), allelic richness (Nar, standardized to a sample of 14), genetic diversity measures (observed heterozygosity [HO] and expected heterozygosity [HE]), inbreeding coefficient (FIS), probability of identity (pI and pIsib), and average relatedness (R; values in boldface type are significantly different from 0, P = 0.001) for contemporary samples of Peary caribou from each region of Ellesmere Island.
When all regions were combined 8 locus pairs showed evidence of significant genotypic disequilibrium. When regions were considered separately fewer loci were statistically linked, and no pair was linked in all regions (NE—RT6/RT9; NW—RT9/BM848; SC—BL42/RT6, RT9/RT24). These loci have not been identified as being linked by previous studies of Peary caribou (Zittlau 2004) or other caribou subspecies (Boulet et al. 2007; Courtois et al. 2003; McLoughlin et al. 2004; Zittlau 2004). These loci were retained for all analyses.
Error rates for microsatellite loci were determined for those sample-locus pairs (e.g., individual 1 at locus RT9) that were profiled in more than 1 polymerase chain reaction (∼20% of total sample loci). The error rate averaged across all loci was 4%, with locus-specific rates ranging from 0% to 8%. Locus-specific error rates were obtained by dividing the number of errors detected (whether due to allelic drop out, false alleles, or scoring error) by the number of replicated sample-locus pairs.
Levels of genetic diversity: historic versus contemporary
Seven haplotyes were identified in the combined data set (Table 2). Two mtDNA haplotypes in the historic NE sample were observed only in 1 specimen each. Because extraction of genetic material from historic NE samples presents a number of technical challenges, we opted to group these 2 haplotypes (haplotypes 6 and 7) with the haplotype that was the most similar and observed in multiple samples (haplotype 5; Table 2). This conservative action resulted in 4 mtDNA haplotypes observed in the historic NE sample, 4 mtDNA haplotypes from the contemporary sample, and 3 mtDNA haplotypes that were shared between the 2 temporal data sets (Table 2). Gene diversity (0.74), mean number of pairwise differences (4.36), and nucleotide diversity (0.02) were highest in the historic NE sample compared to contemporary NE (0.68, 1.24, and 0.006, respectively) or to the entire contemporary data set (0.68, 1.32, and 0.006, respectively).
Based on mtDNA, significant FST differentiation was observed between the historic NE and contemporary data sets due, in part, to the presence of 1 divergent haplogroup identified in the historic NE samples (Table 4). Mitochondrial differentiation was least between historic and contemporary NE samples; however, this differentiation was significant (Table 4, below diagonal). Contemporary NE samples were not significantly differentiated from contemporary CE or SC using mitochondrial markers.
Table 4
Pairwise FST-values among regions computed from mitochondrial control region sequence (214 bp) data. Values below the diagonal are pairwise FST computed from haplotype frequencies only. Values above the diagonal are computed from the matrix of Tamura–Nei (Tamura and Nei 1993) distances. Significant values assessed with 1,000 permutations are indicated in boldface type. Regions of Ellesmere Island are abbreviated as follows; northeast (NE), northwest (NW), central (CE), and south-central (SC).
The haplogroup observed in the historic NE sample was absent from the contemporary sample clusters within Rangifer but away from most Peary caribou haplotypes (results not shown—tree in Appendix III). This could be interpreted as possible contamination. However, given the controls in place to eliminate contamination, the lower discriminatory power due to using a short section of the fragment (214 nucleotides, 11 variable sites), and the well-recognized polyphyly observed in Rangifer when examining larger control region sequences or the cytochrome-b gene (Cronin et al. 2005; Flagstad and Røed 2003; Gravlund et al. 1998), we believe that these represent real haplotypes.
Levels of genetic diversity: contemporary
Mean number of microsatellite alleles per locus for each region ranged from 3.3 in NW to 5.8 in CE (Table 3). Allelic richness estimates indicate a small difference in numbers of alleles when controlled for sample sizes (Table 3). Observed heterozygosity ranged from 0.55 to 0.68 and expected heterozygosity from 0.50 to 0.70 (Table 3). Microsatellite differentiation was low but significant between NE and most other regions (Table 5). We observed no differentiation between contemporary CE and SC (FST = 0.009). Significant differentiation was detected between NW and all other regions (Table 5).
Table 5
Pairwise FST-values calculated using microsatellite data. Boldface values are significant. Regions of Ellesmere Island are abbreviated as follows; northeast (NE), northwest (NW), central (CE), and south-central (SC).
Average relatedness was significantly different from 0 for NW and NE but not for CE or SC. The NW region had significantly higher relatedness compared to the other regions (Table 3). These data showed positive FIS-values in all regions except NW, which had a negative FIS (−0.11). Probability of drawing 2 identical genotypes (pI) from the full data set was in the order of 1 in 3.9 × 106, but this was less for NW (1 in 5.5 × 104). Using the correction for siblings (pIsibs), these odds dropped considerably. In NW this probability of identity was 1 in 138 (Table 3).
Population bottlenecks
The NE region exhibited a signature of a population bottleneck when the infinite alleles model was assumed (P = 0.002), but not under the stepwise mutation model or the 2-phase model (Table 6). A similar result was obtained for CE and the combined data set (Table 6). The NW and SC regions do not show any signature of a bottleneck using the BOTTLENECK program under any mutation model. Cornuet and Luikart (1996) point out that deviations from Hardy–Weinberg equilibrium can affect BOTTLENECK results; when RT24 and RT6 were removed (as per MICRO-CHECKER) the only modification was that NW, under the infinite alleles model, became significant (P = 0.027).
Table 6
Bottleneck results for Peary caribou from each region of Ellesmere Island. M-ratio is reported with interlocus variation in parentheses. Bottleneck results reported as the 1-tailed probability of obtaining the heterozygosity excess observed (P het. excess), based on Wilcoxon's signed-rank test. Values in boldface type are significant at α = 0.05. IAM = infinite alleles model; TPM = 2-phase model; SSM = stepwise mutation model.
The M-ratio for the entire data set was 0.71, which is slightly above the 0.70 threshold indicated by Garza and Williamson (2001) as indicating a recent bottleneck event. When the M-ratio for each region was calculated all regions were below the threshold, with CE having the highest ratio (M-ratio = 0.69) and NW having the lowest ratio (M-ratio = 0.60; Table 6).
Discussion
Evaluation of models
The absence of a historical haplogroup (3 similar haplotypes) from the extensive contemporary sample suggests that the harvest of 233 caribou from NE Ellesmere Island by Robert Peary's expeditions was significant. Although these results should be interpreted in light of the small historic sample and large time interval with no population size data, examination of our data from the large contemporary sample suggests that this severe local reduction could have had a long-lasting impact on the resident animals. Subsequently, caribou dispersed into the NE region from areas to the south. Evidence for this line of reasoning, in addition to the loss of diversity, is the significant differentiation between historic NE and contemporary NE combined with the current similarities between NE and regions to the south but not to the closer NW region. The 100-year interval between historic and contemporary sampling leaves the possibility for many population size fluctuations; however, our data support the recolonization model to a greater extent than the recovery model. Although our data provide only indirect evidence, further genetic recovery could be limited if the diversity observed affects the fitness of these animals. Loss of genetic diversity has been observed in other mammals that were harvested extensively in the 1800s and 1900s (Nyström et al. 2006; Weber et al. 2000). The loss of mitochondrial diversity might reflect a loss of overall genetic diversity that, in turn, can limit the ability of a population to respond to environmental challenges (Frankel and Soulé 1981; Lacy 1997; Lande 1988).
In the NE region a bottleneck signature was identified in the microsatellite data using M-ratio and the infinite alleles model. However, these values should be carefully interpreted because the M-ratio value is close to the empirical threshold and the infinite alleles model is the most likely model to detect a bottleneck. These results suggest there was a bottleneck but do not unequivocally support either of our models, because the bottleneck could be due to recovery by survivors or from a founder event associated with the recolonization model. Immigration, even at low rates, may reduce the signature of a bottleneck for the M-ratio by “filling in the gaps” in the allele frequency range (Garza and Williamson 2001). The transient excess of heterozygotes that BOTTLENECK detects will be lost within a few generations (2Ne to 4Ne generations—Piry et al. 1999).
Although the continued presence of Peary caribou in NE is positive from a conservation perspective, limited recovery suggested by comparing Peary's harvest to current observations, coupled with the observed loss of diversity, should be a warning regarding the trajectory on which overharvest can place species. Direct comparisons of nuclear diversity are not possible, because the markers used differ among studies. However, levels of heterozygosity and mean number of alleles observed in these data generally were lower than observed in caribou from Quebec and Labrador, Canada (Boulet et al. 2007), but higher than observed in a survey of caribou subspecies (Cronin et al. 2003, 2005). Examination of our data provides us with a better understanding of the population genetic structure of Peary caribou on Ellesmere Island and the potential long-lasting impact of a localized harvest event. These types of harvests can be especially significant to large mammal populations that occur at low density and in extreme environments.
Northwest region: isolation
In contrast to NE, the NW region appears to be relatively isolated from adjacent regions. Two unique genotypes (individuals) were obtained from 2 samples collected from the 14 animals seen in 2004, and 17 unique genotypes were obtained from 24 samples collected during a more systematic survey in 2006 when 57 animals were seen. Given the coverage of the 2006 survey of this area (Campbell 2006), it is reasonable to assume that a large portion of the animals were detected. The observation that this sample had 1 predominant mitochondrial haplotype (15 of 17), the highest mean relatedness (0.063), and the lowest observed heterozygosity (0.55), suggests prolonged isolation and at low effective population size. This region also was significantly differentiated from other regions when assessed at microsatellite and mitochondrial markers. These results could arguably be due to the smaller NW sample, yet they differ significantly from the SC region, which had a smaller sample size but is connected to the CE region and to unsampled regions of southernmost Ellesmere Island. The low M-ratio (0.60) suggests that the NW region has been subjected to a severe bottleneck event in the recent past. Garza and Williamson (2001) indicate that the M-ratio will increase if the population size recovers, whereas the number of alleles will continue to decrease for several generations. In contrast, if the population remains low after a bottleneck, the ratio will continue to decline along with the number of alleles. The low allelic diversity, mtDNA haplotype diversity, and M-ratio all suggest a prolonged bottleneck.
If the NW region was completely isolated, we would predict a high degree of inbreeding and a positive FIS. Instead, we observed a negative FIS. Storz (1999) indicated that negative FIS-values are a general characteristic of many mammalian social groups. A harem-polygynous breeding system where few males (usually immigrants) sire most of the young, coupled with a large proportion of offspring making up the sample, will result in elevated levels of heterozygosity (Lawler et al. 2003). Unfortunately, little is known regarding the mating system in Peary caribou. Caribou mating strategies have been described as female-tending in large herds of barren ground caribou (R. t. groenlandicus) to harem-guarding in woodland caribou (R. t. caribou—Banfield 1974; Geist 1998). Observations that would allow us to place Peary caribou into 1 of these categories are currently lacking, and these strategies likely shift in relation to environmental conditions. The negative FIS suggests a harem mating system and large variance in male reproductive success. This could be a result of the relatively isolated area of NW, which may favor male harem-defense or allow females to visit all potential mates and choose the “best.” This mating behavior could have important repercussions for isolated populations similar to those encountered in the artificial isolation that is created by captive breeding and reintroduction programs (Wilson et al. 2005). Negative FIS also could be due to selection for heterozygosity (i.e., balancing selection) on linked functional genes, or inbreeding avoidance. These processes could be occurring, but our data do not allow us to explore these possibilities. Tests of inbreeding avoidance on semidomesticated reindeer (R. t. tarandus) did not find supporting evidence (Holand et al. 2007).
Central and south-central Ellesmere
The weak bottleneck signature (close to the 0.70 threshold) in the CE samples suggests that the severe weather events observed farther south did not affect caribou populations in this area to the same extent. Our M-ratio values are similar to those found in Svalbard reindeer (0.67–0.76—Côté et al. 2002), but are in strong contrast to those reported for 6 islands in the southwestern Canadian Arctic Archipelago (0.35–0.47—Zittlau 2004). Severe winter weather and icing events have been used as a blanket explanation for the dramatic population declines seen in Peary caribou (COSEWIC 2004). Although this is the most parsimonious explanation for the southern islands that were surveyed, it may not be significant across the entire Canadian Arctic Archipelago given the range of environmental regimes from arctic maritime to arctic desert (Kojima 1991). It is possible that the lower overall abundance (Taylor 2005) and greater dispersion of groups (Gauthier 1996; Riewe 1973) on Ellesmere Island could have allowed a greater proportion of the population to survive these events. However, when isolated, a low abundance and high-dispersion system would be susceptible to local events such as overharvest or genetic drift. In either case the extreme declines in the south can serve as warning of the potential for dramatic declines under poor weather conditions.
Conservation implications
Ellesmere Island, with its glaciers and fiords, may naturally fragment Peary caribou groups. Retaining connectivity among fragments will be important to avoid genetic drift or local extinction. The haplotype diversity loss in northern Ellesmere Island between the early 1900s and this contemporary snapshot suggests that the potential to lose genetic variability is real. Since 1900 other areas have experienced more dramatic losses with the extinction of R. t. dawsoni and R. t. groenlandicus (Banfield 1961). Reduced connectivity can have major impacts on genetic diversity. As we suggest for Peary caribou in the NW region, this can be especially significant when isolation is coupled with a mating strategy that further lowers the effective population size. Yet, heterozygosity can sharply increase or be maintained with only an occasional immigrant to the gene pool. Whether it is genetic variability or heterozygosity that is most important for the overall fitness of a population might depend on the threats that population faces. Low genetic variability has been linked to increased extinction rates (Lande and Shannon 1996), whereas heterozygosity has been correlated with various individual fitness traits including parasite resistance (Hansson and Westerberg 2002; Rijks et al. 2008).
Comparisons of genetic diversity between the core and the periphery of a species' range have indicated that diversity generally declines as one approaches the range limit (Garner et al. 2004; Kyle and Strobeck 2002; Vucetich and Waite 2003). However, this trend is not always associated with increased risk of extinction (Channell and Lomolino 2000; Lomolino and Channell 1995), and therefore creating management strategies for populations at the range limit is difficult. In the case of Peary caribou on Ellesmere Island, the range limit is also an absolute barrier, with no suitable caribou habitat beyond the coastline. Therefore, it is difficult to predict the response of this population to future climate changes. As suggested by Moen (2008) in reference to Swedish semidomesticated reindeer, if productivity increases, these areas might support larger numbers of caribou; however, environmental changes also could result in these populations remaining low or declining to a new equilibrium. Continued monitoring will be important to assess these changes because our results for Peary caribou highlight the sensitivity of populations at the range limit to stochastic events and the potential limits to recovery following significant perturbations.
Acknowledgments
This research could not have been accomplished without the hard work of those people and agencies that surveyed Ellesmere Island and collected samples (2004: M. Manseau, T. Mullin, D. Haché, and B. Fyten; 2006: M. Campbell, G. Hope, D. Stern, G. Mouland, and J. Hudson). The 2004 collection was conducted by Parks Canada, with significant support provided by Quttinirpaaq National Park and the Polar Continental Shelf Project. The 2006 collection was a cooperative effort between the Government of Nunavut (Department of Environment), Parks Canada, Agriculture Canada, Polar Continental Shelf Project, the Resolute Bay and Grise Fiord Hunters and Trappers Organizations, and the Nunavut Wildlife Management Board. We are grateful for the assistance of M. Curteanu for help in preparing and shipping the samples and to J. Keeney for assistance in mapping. Thanks also go to the American Museum of Natural History for subsampling and providing hide for DNA extraction. Funding for this project was provided by the Western Canada Service Centre of Parks Canada, the Nunavut Field Unit of Parks Canada, Quttinirpaaq National Park of Canada, and Parks Canada Species at Risk Recovery Action and Education Fund, a program supported by the National Strategy for the Protection of Species at Risk. Additional funding was provided by Natural Sciences and Engineering Research Council to SDP and PJW. Reviews from M. Cronin and C. Krajewski provided helpful comments on the manuscript.
Literature Cited
Appendices
Appendix I
Primer concentrations and locus included in each of the 4 multiplexed reactions. Dropped loci are indicated in boldface type.{ label needed for table-wrap[@id='i1545-1542-91-3-698-t07'] }
Appendix II
American Museum of Natural History (AMNH) specimens and corresponding haplotypes.{ label needed for table-wrap[@id='i1545-1542-91-3-698-t08'] }
Appendix III
Neighbor-joining tree illustrating the polyphyletic nature of a tree derived from a short fragment in a species with a complex evolutionary past. Tree was determined under the Tamura–Nei model (Tamura and Nei 1993) and tested for support with 1,000 bootstrap replicates. Samples include haplotypes identified by Flagstad and Røed (2003), haplotypes observed in this study (•), and unique haplotypes identified from 42 woodland caribou that were part of concurrent studies in the non-paleo laboratories at Trent University and the Natural Resources DNA Profiling and Forensic Centre affiliated with Trent University (possible sources of contamination). Haplotypes are designated by species (Rt = Rangifer tarandus), subspecies (c = caribou, gr = groenlandicus, gt = granti, p = pearyi, and pt = platyrhynchus), GeneBank accession number, and haplotype name.{ label needed for fig[@id='i1545-1542-91-3-698-f301'] }
