Selection and drift influence genetic differentiation of insular Canada lynx (Lynx canadensis) on Newfoundland and Cape Breton Island

Abstract Island populations have long been important for understanding the dynamics and mechanisms of evolution in natural systems. While genetic drift is often strong on islands due to founder events and population bottlenecks, the strength of selection can also be strong enough to counteract the effects of drift. Here, we used several analyses to identify the roles of genetic drift and selection on genetic differentiation and diversity of Canada lynx (Lynx canadensis) across eastern Canada, including the islands of Cape Breton and Newfoundland. Specifically, we assessed whether we could identify a genetic component to the observed morphological differentiation that has been reported across insular and mainland lynx. We used a dinucleotide repeat within the promoter region of a functional gene that has been linked to mammalian body size, insulin‐like growth factor‐1 (IGF‐1). We found high genetic differentiation at neutral molecular markers but convergence of allele frequencies at the IGF‐1 locus. Thus, we showed that while genetic drift has influenced the observed genetic structure of lynx at neutral molecular markers, natural selection has also played a role in the observed patterns of genetic diversity at the IGF‐1 locus of insular lynx.


| INTRODUCTION
Island populations are extreme examples of genetic divergence caused by geographical and landscape barriers to gene flow, where large bodies of water inhibit the immigration of new individuals from the mainland (MacArthur & Wilson, 1967). Such insular populations are often characterized by small population sizes, which are largely the result of founder events and recurrent population bottlenecks, leading to the low-genetic diversity of island-dwelling populations via genetic drift (Frankham, 1997;Nei, Maruyama, & Chakraborty, 1975). Ultimately, small population sizes can result in an increased probability of extinction through demographic stochasticity and inbreeding depression, in addition to limiting the adaptability of island populations to environmental change (Hedrick & Kalinowski, 2000;Reed & Frankham, 2003).
For example, a decrease in genetic diversity due to genetic drift could impede adaptation to new environmental conditions by reducing standing adaptive genetic variability (Agashe, Falk, & Bolnick, 2011;Bridle & Vines, 2007;Swaegers et al., 2013). Additionally, insular populations are likely to be prevented from adapting by means of evolutionary rescue, in which immigrants promote adaptive evolution by introducing novel, beneficial alleles lost to drift, increasing population size and thus mutational opportunities, and reversing negative growth rates (Bell, 2013;Bell & Gonzalez, 2011).
Apart from their generally high susceptibility to extinction, island populations have long been a focus of study in the evolutionary literature due to their often unique contribution to the adaptive diversity of a species (Watanabe, Kazama, Omura, & Monaghan, 2014). For example, some characteristics appear to evolve faster than expected on islands because high standing adaptive genetic variation is believed to be the most influential mechanism for adaptation (Barrett & Schluter, 2008), and island populations generally exhibit low levels of genetic diversity (Nei et al., 1975). Such fast-paced evolution on islands is likely determined by strong selective pressures for specific traits, driven by the highly divergent environments and community structures often encountered on island landscapes that are necessary to overwhelm the effects of genetic drift. This adaptive divergence is also likely reinforced by a lack of homogenizing gene flow with mainland populations. For example, Millien (2006) compared the rates of evolution between island and mainland populations of a number of mammal species and found that significant morphological changes can occur over short time intervals on islands. These changes are most commonly observed as differences in size (Lomolino, 1985(Lomolino, , 2005Lomolino et al., 2012Lomolino et al., , 2013, where smaller mammals on islands evolve toward gigantism (e.g., Millien & Damuth, 2004) and larger mammals evolve toward dwarfism (e.g., Roth, 1992;Sondaar, 1991). This graded trend, both across and within island taxa, was termed the "island rule" by Van Valen (1973) and has since been supported by numerous studies comparing species inhabiting both islands and the mainland (e.g., Dayan & Simberloff, 1998;Heaney, 1978;Kranzowski, 1967). The generality of the island rule is driven by the effect of body size on immigration potential, resource requirements, and inter-and intraspecific interactions of island species (Lomolino, 2005). Larger species, for example, experience a limitation in the amount of resources on islands and thus increased intraspecific competition, resulting in the higher fitness of smaller individuals who require less energy for survival and reproduction (Case, 1978;Lomolino, 1985Lomolino, , 2005. The validity of the island rule in regard to carnivores, however, has been challenged by some (i.e., Meiri, Dayan, & Simberloff, 2004, 2006, and the molecular mechanisms that underlie morphological differentiation are largely unknown. As ecologically marginal and unique populations (e.g., edge or insular populations) are likely to encompass important components of genetic diversity, the loss of such populations are likely to result in the loss of rare and potentially beneficial adaptive genotypes (Watanabe et al., 2014). Such implications extend to management scenarios where human intervention is used as a means for the evolutionary rescue of island populations. If unique adaptations are present on islands, outbreeding depression may result from the introduction of maladapted individuals; the influx of their genes into the island population may swamp the well-adapted gene pool and break down coadapted gene complexes (Miller, Poissant, Hogg, & Coltman, 2012).
On the other hand, if genetic drift is predominating the distribution of island alleles and genotypes, the introduction of outside individuals may be an effective management strategy. Wide-ranging species present unique challenges in this regard, due to high rates of gene flow and resulting genetic panmixia of the species across their range.
For such species, insular populations may thus represent the only reservoirs of unique genetic material due to their isolation from the mainland, occupation of unique ecological settings, and exposure to divergent selection pressures (Gaston, 2003).
A classic example of a wide-ranging and panmictic species is the Canada lynx (Lynx canadensis), which occurs throughout the boreal forests of Canada and Alaska, USA, and parts of the lower 48 United States. Canada lynx also inhabit the islands of Newfoundland and Cape Breton Island, Canada. Little is known, however, about the origin and subsequent isolation of these populations in comparison with the mainland population (Khidas, Duhaime, & Huynh, 2013). Row et al. (2012) conducted a genetic survey of lynx across their Canadian range, and found that lynx in Newfoundland were highly distinct at 21 neutral microsatellite markers. More recently, Koen, Bowman, and Wilson (2015) identified genetic differentiation between Canada lynx located north and south of the St. Lawrence River in Quebec, Canada. No such survey has been conducted for lynx on Cape Breton Island.
In addition to the inherent isolation of lynx on Cape Breton Island, morphological differentiation of Canada lynx across its range is greatest between individuals on Cape Breton Island and mainland Canada (Khidas et al., 2013). Specifically, of 18 craniodental traits measured by Khidas et al. (2013), both male and female lynx from Cape Breton Island had smaller and more variation in skull shape compared with the mainland population. This suggests that the body size of lynx on islands is smaller than in the mainland population, although this has not been explicitly tested. If consistent with the island rule, the smaller size of lynx on Cape Breton Island could be the result of their historically high population densities on the island (10-20 lynx/100 km 2 during the peak of the population cycle; Parker, 1981Parker, , 2001, resulting in increased intraspecific competition for limited resources. A similar negative relationship between population density and body size has been reported for lynx in Alaska (Yom-Tov et al., 2007). Further, evidence of differences in the social behavior (lynx on Cape Breton Island have been reported to regularly socialize in spruce bogs; Parker, 1981) and prey consumption (lynx on Cape Breton Island consume higher proportions of white-tailed deer carrion or bait (Parker, Maxwell, Morton, & Smith, 1983) versus mainland lynx which rely on red squirrel as a secondary food source ) supports evolutionary differentiation between Cape Breton and mainland lynx populations.
Given the preceding evidence that lynx on Cape Breton Island have diverged phenotypically from individuals on the mainland, it is reasonable that both genetic drift and selection have contributed to their evolution and differentiation since their establishment on the Island, and perhaps on Newfoundland as well. Here, we used both neutral microsatellites and one non-neutral (functional) genetic marker (insulin-like growth factor-1; IGF-1) to evaluate the influence of genetic drift and natural selection on insular and mainland populations of Canada lynx in eastern Canada. First, we were interested in determining whether lynx on Cape Breton Island are genetically differentiated from lynx on mainland Canada. Genetic distinctness is likely for this population, given: (1) genetic differentiation at neutral molecular markers is proportional to time since divergence under an island model when it can be assumed that the population has descended from an ancestral population and since diverged without gene flow (Nielsen & Wakeley, 2001), and (2) lynx are thought to have occurred on Cape Breton Island as early as the 1600s (although these reports refer to Nova Scotia and not Cape Breton Island specifically; Denys, 1672). Second, as the function of the IGF-1 gene has been linked to body size in mammals (e.g., Sutter et al., 2007), we were interested in determining whether differentiation may be observable at the CA dinucleotide repeat between island and mainland lynx groups for which known morphological differentiation exists (Khidas et al., 2013).
This repeat lies in the promoter region of the IGF-1 gene, containing regulatory elements that may alter transcription of the gene (Wagner et al., 2005). After Khidas et al. (2013), we hypothesized that lynx on Cape Breton Island are genetically divergent from mainland lynx, and predicted that we would observe high genetic distinction between lynx on Cape Breton Island and on the mainland at both neutral microsatellite markers and the non-neutral IGF-1 dinucleotide repeat. Third, we were interested in assessing the contributions of genetic drift and natural selection to the observed allele frequencies of the IGF-1 locus in island versus mainland lynx. Lastly, we wanted to test whether genotypes at the IGF-1 locus were associated with measures of morphology in lynx on Cape Breton Island. If the gene is under selection or linked to a nearby locus under selection in lynx, we predicted that we would observe a relationship between common IGF-1 alleles found in lynx on Cape Breton Island and morphology.

| Study area
The mainland Nova Scotia population of Canada lynx was extirpated in the 1920s-1930s (Smith, 1940;de Vos & Matel, 1952), and now, the only remaining population of lynx in Nova Scotia occurs in the western plateau of Victoria and Inverness counties in the highlands of Cape Breton Island (Quinn & Parker, 1987; Figure 1). Cape Breton Island has a land area of approximately 9,500 km 2 and is separated from mainland Nova Scotia by the Strait of Canso, which is approximately 1.6 km wide at its narrowest point (Matrin, 1837; Figure 1).
Historically, population densities of lynx on Cape Breton Island were high (Parker et al., 1983); however, more recent estimates of lynx density and area of occupancy found that only 50-500 animals remained on the island, depending on the stage of the population cycle (Parker, 2001 Uncertainty remains, however, as to whether migration to and from the island results in successful genetic exchange with mainland lynx (Poole, 2003).

| Sample collection and DNA extraction
We sampled Canada lynx from Quebec (both north and south of the St. Lawrence River), New Brunswick, Labrador, and the islands of Newfoundland and Cape Breton (Figure 1). With the exception of Cape Breton Island, the samples we used in this study were the same as those in Koen et al. (2015), with some small changes in sample size.
We collected hide samples from lynx that were harvested for their fur
We further amplified the dinucleotide microsatellite repeat within the IGF-1 functional gene in all samples with the primers described by Kirkpatrick (1992). For the IGF-1 locus, we conducted amplification in a 10 μl reaction containing deionized water (Invitrogen), 1X PCR Reaction Buffer (Invitrogen), 2 mmol/L MgCl 2 (Invitrogen), 0.2 mmol/L dNTP solution (Invitrogen), 0.2 mg/ml BSA, 0.4 μmol/L forward and reverse primers (forward primer labeled with the fluorescent dye HEX) (Integrated DNA Technologies), 0.025U Taq DNA Polymerase (Invitrogen), and 5 ng of DNA. We ran the PCR in a Bio-Rad DNA Engine Dyad Disciple thermocycler under the following conditions: 95°C for 10 min followed by 30 cycles of 94°C for 30 s, 62°C for 1 min, and 72°C for 1 min, and completed with a step of 65°C for 15 min.
For genotyping of the IGF-1 locus, we mixed 5 μl of size standard (MapMarker 1000 X-Rhodamine, BioVentures) into 1 ml of HiDi Formamide (Applied Biosystems), and we added 9.5 μl of this product to 0.5 μl of each amplified sample.
We performed genotyping on the Applied Biosystems 3730 DNA Analyzer. We scored genotypes for neutral microsatellites in lynx from Cape Breton Island and all lynx at the IGF-1 locus with Softgenetics LLC GeneMarker AFLP/Genotyping Software Version 1.91 (State College, PA, USA). We used manual scores from Koen et al. (2015) for all other samples. We omitted samples that were missing alleles at more than 2 of the 14 microsatellite loci. We analyzed the IGF-1 locus separately. We used the R (R Core Team 2016) package adegenet, version 2.0.1 (Jombart, 2008)

to estimate expected (H E ) and observed (H O )
heterozygosity for each locus, and conducted a t-test to evaluate whether average expected and observed heterozygosity were significantly different across loci. We also used this package to conduct Bonferroni-corrected (α = .05) chi-squared tests to detect loci that departed from Hardy-Weinberg equilibrium (HWE). We used the software GENEPOP (version 4.2; Raymond & Rousset, 1995;Rousset, 2008) to conduct Bonferroni corrected (α = .05) tests for linkage disequilibrium (LD). We used HP-RARE version 1.1. (Kalinowski, 2005) to estimate allelic richness (AR) for each of our 6 sites, and corrected our estimates with rarefaction for our lowest sample size (N = 14). We used the R package diversity (Keenan, McGinnity, Cross, Crozier, & Prodöhl, 2013) to obtain multiple estimates of genetic differentiation (F ST [Weir & Cockerham, 1984]; G ST [Nei, 1973]; and Jost's D [Jost, 2008]) as well as F IS , with 95% confidence intervals (CIs) on each estimate (999 bootstraps). We estimated AR, genetic differentiation, and Breton Island (N = 38). All of these samples were from the same individuals that we sampled for neutral markers.
We used the Bayesian-clustering software STRUCTURE version 2.3.4 (Pritchard, Stephens, & Donnelly, 2000) to identify genetic clusters based on our 14 neutral microsatellite loci (excluding the IGF-1 locus). We tested scenarios of K = 1-9 and ran 10 repetitions for each value of K. We used a burn-in length of 500,000 Markov chain Monte Carlo iterations, followed by a run length of 1,000,000 iterations.
We selected an admixture model without prior location information and used a correlated allele frequency model. We selected the most likely number of genetic clusters using the Evanno method (Evanno, Regnaut, & Goudet, 2005) implemented in the software STRUCTURE  (Earl & vonHoldt, 2012). We summarized our 10 replicate runs for each value of K with the software CLUMPP (Jakobsson & Rosenberg, 2007), and visualized them with the software DISTRUCT (Rosenberg, 2004). To add support to our STRUCTURE results, we also estimated pairwise measures of F ST , G ST , and Jost's D between clusters identified in STRUCTURE, and implemented Fisher's exact tests to estimate whether genotype distributions at our 14 neutral microsatellites were significantly different between each pair of clusters (2000 MC replications, p < .001) with the R package diversity (Keenan et al., 2013). We assigned admixed individuals to the cluster that accounted for more than 50% of their ancestry.

| Assessing the influence of genetic drift versus selection at the functional IGF-1 locus
There are two alternative hypotheses explaining genotypic and allelic variation in insular lynx populations. First, allelic distributions may be largely influenced by genetic drift. If this were the case, we would expect to find evidence of small effective population sizes (N e ), a history of population bottlenecks, and/or reduced genetic diversity across loci in island but not mainland populations. Alternatively, insular populations could be experiencing strong directional or positive selection, where few alleles are driven to high frequencies on the islands (not necessarily the same alleles as those at high frequency on the mainland) due to differences in environmental conditions, prey base and/or inter and intraspecific competition. In this case, we would expect values of genetic differentiation that are greater at the loci under selection (e.g., IGF-1) in comparison with neutral markers (e.g., F ST outliers), or a correlation of alleles with environmental or other variables that differ between the examined mainland and insular populations (e.g., morphology). To differentiate between these hypotheses, we conducted multiple analyzes.
To identify the contribution of genetic drift, we followed the methodology of Funk et al. (2016) using our neutral microsatellite marker dataset. First, to describe genetic diversity, we calculated three measures of within-population genetic variation; observed (H O ) and expected (H E ) heterozygosity, and allelic richness (AR). We estimated these metrics using the same approaches as for our regionally defined populations (AR was corrected for our lowest sample size [N = 26] using rarefaction). We also estimated N e for each of our genetic clusters using the molecular co-ancestry method of Nomura (2008), as implemented in NeEstimator (version 2.01; Do et al., 2014).
We used this method because it has been shown to provide unbiased estimates of N e , estimated as the effective number of breeders (N eb ) for noninbred populations, and its narrower confidence intervals are more practical for interpretation (Nomura, 2008). Further, we estimated whether any of our genetic clusters have experienced recent bottlenecks using the software BOTTLENECK (version 1.2.02; Cornuet & Luikart, 1996). We used the default parameters under the two-phase model (TPM) of mutation for our loci, as microsatellites with <3 bp repeats (e.g., our neutral dinucleotide repeats) generally have mutation models that include multiple-step mutation events (Di Rienzo et al., 1994), and are therefore somewhere in between the standard step-wise mutation model (SMM) and the infinite alleles model (IAM; Cornuet & Luikart, 1996). We ran 1,000 iterations for each test and used one-tailed Wilcoxon signed rank tests to detect significant excess heterozygosity that is generally indicative of a bottleneck. Lastly, if genetic drift was influencing insular lynx, we would expect pairwise differentiation between genetic clusters to be negatively correlated with estimates of within-population genetic variation (used as an index of historical genetic drift). Such an analysis has been used to identify genetic drift in populations of reptiles (Jordan & Snell, 2008), fish (Whiteley et al., 2010), and mammals (Funk et al., 2016). To test this hypothesis, we conducted linear regressions between pairwise estimates of genetic differentiation (F ST and Jost's D) within our 4 genetic clusters, and average measures of withinpopulation genetic diversity (H O , H E , and AR) at our 14 neutral microsatellite markers. As gene flow with mainland populations is unlikely to have significantly contributed to within-population genetic diversity of the islands (Funk et al., 2016), estimates of within-population genetic diversity are an ideal index of historical genetic drift for this analysis.
We used the IGF-1 locus as an index of selection in our insular lynx populations. We plotted allele frequencies for our four genetic clusters (as identified with neutral markers by the program STRUCTURE; see section 3) at the IGF-1 locus to visualize the distribution of alleles, and used Fisher's exact tests implemented in GENEPOP (version 4.2; Raymond & Rousset, 1995;Rousset, 2008) to estimate whether allelic distributions were significantly different between each cluster. We also used a coalescent-based (Beaumont & Nichols, 1996) outlier detection approach implemented in the software LOSITAN (Antao, Lopes, Lopes, Beja-Pereora, & Luikart, 2008), which uses F ST as a function of heterozygosity to identify loci under selection. We calculated the "neutral" mean F ST , where we first simulated removal of potentially selected loci prior to computing the initial mean F ST upon which putative adaptive loci are identified. We also selected the option to "force mean F ST ," in which we approximated a more precise F ST by running a bisection over repeated simulations. We ran 50,000 simulations at a 95% confidence interval, and selected a stepwise mutation model, which is commonly used to describe tandem repeat markers (Antao et al., 2008;Fan & Chu, 2007). Our LOSITAN tests were conducted using the genetic clusters identified with STRUCTURE. We ran an initial test across all four clusters, and then conducted tests removing each genetic cluster iteratively with replacement. In cases where LOSITAN identified a signature of selection at any locus, two additional independent tests were conducted as confirmation. We also plotted the pairwise Jost's D values of our neutral marker dataset with standard error of the mean intervals, and then overlaid each pairwise estimate of Jost's D at IGF-1 to determine whether these estimates fell within or outside of our expected range of neutrality. We used Jost's D here, as a complementary analysis to LOSITAN because F ST and its analogs (e.g., G ST ) are dependent on, and can be biased by within-population heterozygosity, especially in cases when migration rates are expected to be low (Hedrick, 1999), whereas Jost's D uses the effective number of alleles (Kimura & Crow, 1964) which is a more appropriate measure for our outlier-based approach.

| Testing for associations between IGF-1 alleles & morphology
As significant morphological divergence has been documented in insular lynx in comparison with mainland populations (Khidas et al., 2013), we were interested in testing whether these differences in morphology might be observable at a microsatellite repeat within the IGF-1 gene that has been linked to body size in some mammals (e.g., Baker, Liu, Robertson, & Efstratiadis, 1993;Sutter et al., 2007) We obtained seven morphological measures for each of our individual lynx genotyped on Cape Breton Island: total length, length of tail, length of hind foot, weight of skinned carcass, total length of skull, zygomatic width, and mandible length. We conducted 2-way ANOVAs on each of these measured traits to determine whether mean trait values differed significantly by IGF-1 genotype, and corrected p-values using Bonferroni correction (α = .05). Individuals with missing measurements were removed. We controlled for sex by including it as our second factor in each analysis, but were unable to control for age due to sample size limitations. Thus, to indirectly assess the influence of age on our measured traits, we estimated correlations between age and each trait using Pearson's correlation coefficients (for normally distributed data), or Spearman's rank correlation coefficients (for nonnormally distributed data).

| Genetic diversity and population structure
The number of alleles per locus ranged between 5 and 21, and mean expected and observed heterozygosity were significantly different (p < .001) across our set of 15 neutral and functional loci, suggesting deviations from HWE (Table S1) River) that significantly indicated inbreeding (95% CIs did not overlap 0), but the range of values of F IS was small (0.005-0.038 and 0.023-0.079 for Quebec north and south of the St. Lawrence River, respectively), and thus inbreeding is likely not prevalent in these populations (Table S2)

. Allelic richness was lower in Newfoundland and Cape
Breton Island in comparison with all mainland sites at both neutral markers (Newfoundland AR = 3.52; Cape Breton Island AR = 2.64) and the IGF-1 locus (Newfoundland AR = 2.00; Cape Breton Island AR = 2.00; Table S2). Private allelic richness, however, was lower in sites south of the St. Lawrence River (i.e., Quebec south of the St. Lawrence River, New Brunswick, and Cape Breton Island) than sites north of the river (i.e., Quebec north of the St. Lawrence River, Labrador, and Newfoundland) at both marker sets, with the exception that Newfoundland had a private allelic richness consistent with sites south of the St. Lawrence River at the IGF-1 locus ( Table S2).
Estimates of genetic differentiation showed largely similar patterns for each of the metrics we measured (F ST , G ST and Jost's D). Following Keenan et al. (2013), we plotted the relationship between mean number of alleles and each of our genetic differentiation metrics. The positive relationship that we observed between the mean number of alleles and both F ST and G ST indicated that high mutation rates were not biasing these estimates. Thus, we report genetic differentiation as F ST . As we predicted, the most highly differentiated sites were the islands ( Table 1). The next most differentiated pairs of sites involved one of the two insular populations and one of the mainland sites, and differentiation was greater for comparisons including Newfoundland than Cape Breton Island (Table 1) Table 1). There was a slight change in the genetic structure of our sampled sites when genetic differentiation was measured at the IGF-1 locus ( Table 1).
The two insular populations were still the most differentiated sites; however, pairwise comparisons with Cape Breton Island were generally more highly differentiated at this locus than comparisons with Newfoundland. Evidence for higher genetic differentiation of mainland sites on opposite sides of the St. Lawrence River than sites on the same side of the river was also present at this locus. Both Evanno's Delta K and rate of change plots suggested two genetic clusters using our neutral microsatellite data set ( Fig. S1A-C (Table S3). We also found that the distribution of genotypes across our four genetic clusters were significantly different (p < .001) supporting four versus two genetic clusters. In concert with published evidence that Newfoundland (Row et al., 2012) and the St. Lawrence River (Koen et al., 2015) are impediments to gene flow for Canada lynx, we suggest that four genetic clusters are most likely for these data (Figure 1). We found that Newfoundland had greater differentiation at neutral markers and Cape Breton Island had greater differentiation at the IGF-1 locus whether we estimated genetic differentiation across our six sampling sites or four genetic clusters.

| Assessing the influence of genetic drift versus selection at the IGF-1 locus
Estimates of genetic diversity within each of the four genetic clusters were lower for island clusters compared to mainland clusters at our 14 neutral microsatellite markers, as expected under the hypothesis of genetic drift (Table 2). Effective population sizes (N e ) were small for the insular populations, whereas we estimated N e of mainland sites as infinite, and thus likely to be large (Table 2). We detected evidence of bottlenecks in the form of heterozygosity excess for all four genetic clusters (p < .05), while there was slightly higher evidence for bottlenecks in the two mainland clusters (p < .01) than the island clusters (p < .05), there was also greater power in our analyses of the mainland sites as they had much larger sample sizes overall. Linear regressions between pairwise genetic differentiation and each of our within-population measures of genetic diversity (H O , H E , and AR) were T A B L E 1 Pair-wise F ST (Weir & Cockerham, 1984)

| Testing for associations between IGF-1 alleles and morphology
Our Fisher's exact tests indicated that the distribution of genotypes at the IGF-1 locus was significantly different across all pairs of sites (p < .05) excluding some comparisons involving two mainland sites (Labrador and New Brunswick, New Brunswick and Quebec north, and New Brunswick and Quebec south). When we conducted this estimation across genetic clusters identified using STRUCTURE, genotypic distributions were significantly different across all pairwise comparisons (p < .01).
All seven datasets used in our ANOVA analyses, had homogeneous variance. Sex was found to be a significant factor in two of our ANOVA analyses (total length, F = 14.75, df = 1, p = .0005; length of hind foot, F = 11.12, df = 1, p = .0022); however, the factor "IGF-1 genotype" was not significant for any of our traits following Bonferroni correction (Table S4). Age and trait variables were correlated for one data set only (zygomatic width, t = 2.508, n = 34, p = .0176).

| DISCUSSION
We found that lynx on Cape Breton Island were genetically distinct, supporting the hypothesis that the Strait of Canso is a barrier to gene flow for this insular population. Supported by previously published evidence of the Strait of Belle Isle isolating Newfoundland lynx (Row et al., 2012), and the St. Lawrence River isolating lynx north and south of the river (Koen et al., 2015), water barriers seem to be the most sig- Both genetic drift and natural selection seem to influence the genetic differentiation and diversity of insular lynx. Analyses with our neutral molecular markers showed evidence of (1) lower genetic T A B L E 2 Observed heterozygosity (H O ), expected heterozygosity (H E ), and allelic richness (AR) estimates calculated for 14 neutral microsatellite markers and the IGF-1 locus, and effective population size estimates (N e ) with 95% confidence intervals (calculated from 14 neutral microsatellite loci) for 591 Canada lynx (Lynx canadensis) sampled from eastern Canada Lawrence River than in Newfoundland (Environment Canada, 2016).

Genetic cluster H O (neutral) H O (IGF-1) H E (neutral) H E (IGF-1) AR (neutral) AR (IGF-1)
We detected evidence of bottlenecks in our mainland genetic clusters at a higher level of significance than our island groups. We were expecting higher power in our analysis of the mainland groups; however, because the population sizes of our islands are small and thus bottlenecks are not likely to be as strong, reducing the period during which weak bottlenecks can be detected with high power (Cornuet & Luikart, 1996). In general, our BOTTLENECK analysis likely lacked power because the approach has lower power when loci evolve under a SMM (as opposed to the IAM) and have high heterozygosity (0.3 < H < 0.8) (Cornuet & Luikart, 1996). The bottlenecks detected in our mainland groups may merely reflect the difference in cyclic amplitude of Canada lynx (O'Donoghue et al., 1998;Stenseth et al., 1998;Yan, Stenseth, Krebs, & Zhang, 2013), which is lower or absent in peripheral populations of lynx (e.g., edge or island populations) in comparison with core populations (Murray, Steury, & Roth, 2008;Poole, 1994).
In addition to genetic drift, we were able to identify a functional locus that has likely been influenced by natural selection in insular lynx, IGF-1, based on: (1) allele frequency differences; (2)  (3) outlier tests using Jost's D. Combined, these results suggest that selection has also contributed to the genetic differentiation and diversity observed in island lynx populations. While it is true that the lower genetic diversity (e.g., allelic richness) of the IGF-1 locus on island compared to mainland groups of lynx could be the result of genetic drift, selection for fewer, but more favorable alleles on islands are better supported by our cumulative data. This is because both islands have retained the same two alleles despite likely being founded by genetically divergent populations on the mainland (north versus south of the St. Lawrence River). Further, the more common of the two alleles on Cape Breton Island, which is also found at a relatively high frequency in Newfoundland, was not found at high frequencies in either of the mainland sites of eastern Canada that we sampled, and there is likely a lack of gene flow between the two islands (supported by high estimates of F ST and Jost's D). Thus, the unlikely prospect that an uncommon allele on the mainland would be present in the founding populations of both Cape Breton and Newfoundland and subsequently drifted to high frequencies on both islands lends support to the hypothesis that convergent selection has caused the increase in frequency of this allele in both populations of insular lynx simultaneously.
Patterns of differentiation also support a selection-driven mechanism for the observed IGF-1 genotypes in our insular lynx populations. Specifically, while Newfoundland was more highly differentiated than Cape Breton Island with neutral genetic markers, the opposite trend was observed at the IGF-1 locus, suggesting that different forces are influencing genetic differentiation at these markers (i.e., drift vs. selection). Greater neutral genetic differentiation of Newfoundland is consistent with the premise that gene flow with the mainland is less prevalent in Newfoundland versus Cape Breton Island, which exchanges more migrants with the mainland via the Canso causeway (Nova Scotia Lynx Recovery Team 2006). At IGF-1, however, Cape Breton Island was more highly differentiated, which could be attributed to selection occurring over a longer interval time on Cape Breton, as lynx are thought to have inhabited Cape Breton Island since at least the seventeenth century (Denys, 1672), whereas they were first reported on Newfoundland in 1897 (Bangs, 1897). Alternatively, there may be stronger selection on Cape Breton due to increased competition and predation (Parker, 2001;Parker et al., 1983;Vashon, Vashon, & Crowley, 2003), thereby limiting resources and selecting for the smaller body sizes observed on Cape Breton Island. While we were unable to find evidence that the IGF-1 locus was under selection using the F ST outlier detection approach, we were able to identify Cape Breton Island and Newfoundland as outliers using our Jost's D method in place of F ST . This is likely because Jost's D has been shown to perform better at estimating the differentiation in allele frequencies among populations, whereas F ST is preferential for describing the effect of demography on genetic variation (Meirmans & Hedrick, 2011). In our analysis, Cape Breton and mainland group comparisons were significantly more differentiated at the IGF-1 locus than neutral expectations, indicating that divergent selective pressures are driving the differentiation of lynx on Cape Breton Island relative to mainland lynx at this locus. In contrast, lynx from Newfoundland compared to lynx from Cape Breton and Quebec south of the St. Lawrence River/ New Brunswick were significantly less differentiated at the IGF-1 locus than expected from our neutral markers, indicating that similar selective pressures are decreasing differentiation between Newfoundland and other peripheral and isolated populations of lynx. Importantly, genetic differentiation of the mainland groups was very similar and nonsignificant between neutral markers and IGF-1, suggesting that gene flow may overwhelm selection, and that selection on this locus is not strong enough to counteract this effect in these populations, resulting in patterns of differentiation at the IGF-1 locus that mimic neutral processes. Lastly, a nonsignificant deviation between IGF-1 and neutral marker estimates of Jost's D for Newfoundland, and its closest mainland relative, Quebec north of the St. Lawrence River/Labrador, suggests a strong effect of gene flow counteracting selection between these populations. This is more difficult to interpret; however, as gene flow is likely very limited between these two populations, which are separated by the Strait of Belle Isle, a minimum of 17.6 km wide barrier ( Figure 1).
Although we did not find associations between morphological trait measurements and IGF-1 genotype of lynx on Cape Breton Island, there were limitations in our approach to detecting associations between genetic and morphological data. First, we cannot reject the possibility that unmeasured traits may be significantly associated with the IGF-1 genotypes observed on Cape Breton Island. Further, associations between morphology and the IGF-1 locus could be better understood by comparing across groups (e.g., on both islands as well as mainland groups), and increasing sample sizes within groups.
Indeed, the probability of detecting significant associations between morphological traits and IGF-1 alleles would be greater if a larger variation of morphologically distinct phenotypes were compared across island and mainland sites; however, we did not have access to morphological measurements for the mainland lynx samples used in our study.
Coupled with the fact that the Cape Breton Island population has low levels of variation in morphological trait measurements and high levels of gene flow, our analysis likely lacked power to associate common IGF-1 alleles with phenotype. Ultimately, addressing these limitations may allow for a better understanding of the genetic basis of the island rule in mammals. Taken together, our analyses of the IGF-1 locus indicated that while genetic drift may be strong in island populations, the strength of selection may also be strong enough to allow for the adaptive divergence of important traits, including body size.
Disentangling the contributions of genetic drift and natural selection in island populations is an important conservation goal with relevance to species management. The low genetic diversity of island Canada lynx populations via genetic drift reported here may suggest that lynx are limited in their capabilities to track changes in the environment, especially at a time-scale that is consistent with climate change. Specifically, the ability of insular populations to adapt to environmental stressors is largely dependent on the population's standing adaptive genetic variation (Barrett & Schluter, 2008), which is essentially hindered in insular lynx due to their small effective population sizes and overall low levels of genetic diversity. This suggests that these island populations may rely largely on de novo mutations for adaptation, a mechanism that may not be very effective in short time frames. Consequently, the small, insular lynx populations of Newfoundland and Cape Breton Island are likely at higher risk of extirpation (see Frankham, 1997), and thus warrant appropriate conservation efforts.
Although extinction risk is more likely driven by genetic drift than selection in small populations (Whitlock, 2000), the identification of adaptively divergent traits further suggests that insular populations of Canada lynx warrant conservation efforts to preserve the unique alleles that are present at high frequencies on the islands but are largely absent from the mainland. Knowledge of adaptive genetic diversity may also benefit management plans by informing the potential for genetic rescue of these populations. A greater understanding of a larger breadth of adaptive genetic differences and/or similarities, seek to examine a wider variety of adaptive traits to more fully understand the adaptive landscape of both island and mainland lynx.

| CONCLUSIONS
The protection of island populations is often a conservation goal due to their unique genetic contribution to the species. Integral to this goal is understanding the relative roles of genetic drift and natural selection in shaping the genetic diversity and differentiation of insular populations. The protection of adaptive genetic diversity within species is critical for ensuring species persistence in the face of both short-term environmental stressors and long-term climate change.
Additionally, understanding the impact of genetic drift in isolated and often small island populations is a critical step toward understanding whether the loss of genetic variation imposes a significant threat to the long-term persistence of populations via inbreeding depression.
We showed that while genetic drift has played a role in shaping the observed genetic structure of insular Canada lynx populations, adaptive differentiation also seems to contribute to their unique gene pool.
Extending the analysis of genes putatively under selection (Kirk & Freeland, 2011) will allow for a better understanding of the influence of selection on islands in traits that are not phenotypically apparent and will greatly contribute to our understanding of the genetic basis for selection in island environments. This could also translate to a better understanding of the effects of landscape fragmentation on the mainland, which essentially creates "island" populations in species with limited dispersal abilities.