Killer whales (Orcinus orca) in Iceland show weak genetic structure among diverse isotopic signatures and observed movement patterns

Abstract Local adaption through ecological niche specialization can lead to genetic structure between and within populations. In the Northeast Pacific, killer whales (Orcinus orca) of the same population have uniform specialized diets that are non‐overlapping with other sympatric, genetically divergent, and socially isolated killer whale ecotypes. However, killer whales in Iceland show intrapopulation variation of isotopic niches and observed movement patterns: some individuals appear to specialize on herring and follow it year‐round while others feed upon herring only seasonally or opportunistically. We investigated genetic differentiation among Icelandic killer whales with different isotopic signatures and observed movement patterns. This information is key for management and conservation purposes but also for better understanding how niche specialization drives genetic differentiation. Photo‐identified individuals (N = 61) were genotyped for 22 microsatellites and a 611 bp portion of the mitochondrial DNA (mtDNA) control region. Photo‐identification of individuals allowed linkage of genetic data to existing data on individual isotopic niche, observed movement patterns, and social associations. Population subdivision into three genetic units was supported by a discriminant analysis of principal components (DAPC). Genetic clustering corresponded to the distribution of isotopic signatures, mtDNA haplotypes, and observed movement patterns, but genetic units were not socially segregated. Genetic differentiation was weak (F ST < 0.1), suggesting ongoing gene flow or recent separation of the genetic units. Our results show that killer whales in Iceland are not as genetically differentiated, ecologically discrete, or socially isolated as the Northeast Pacific prey‐specialized killer whales. If any process of ecological divergence and niche specialization is taking place among killer whales in Iceland, it is likely at a very early stage and has not led to the patterns observed in the Northeast Pacific.

Understanding patterns of population genetic structure and the processes driving those patterns is highly relevant for conservation purposes, since smaller population units are more vulnerable to extinction (Stevick et al., 2006).
Contrarily to the Northeast Pacific, although dietary variation and some degree of ecological divergence have been reported among North Atlantic killer whales (Foote et al., 2013(Foote et al., , 2011Foote, Newton, Piertney, Willerslev, & Gilbert, 2009;Foote, Vester, Víkingsson, & Newton, 2012), the link between genetic divergence and resource specialization is less clear. In Iceland, not all killer whales appear to specialize on herring and follow it year-round.
Photo-identification of killer whales in Icelandic herring overwintering-and summer-spawning grounds showed that some individuals are sighted in both summer and winter seasons but others are seen only seasonally . Also, some individuals observed only in the winter season were seen moving to Scotland in the summer, where they were seen feeding upon marine mammals (Samarra & Foote, 2015). Isotopic analyses of biopsy sampled individuals with different observed movement patterns show that these largely correspond to different isotopic niche widths (Samarra, Vighi, Aguilar, & Víkingsson, 2017). Individuals seen in both seasons had lower nitrogen stable isotope ratios ( 15 N/ 14 N, represented as δ 15 N), consistent with a diet predominantly composed of herring, while killer whales only seen seasonally (including sampled individuals that travel between Iceland and Scotland) exhibited larger variation in δ 15 N, suggesting that some individuals have a diet including other prey (Samarra, Vighi, et al., 2017). However, there is no social isolation between individuals with different observed movement patterns and isotopic signatures, that is, putative herring-specialists remaining year-round in Iceland have been photographed in close proximity with Icelandic-Scottish killer whales . It is unknown whether the apparent absence of social isolation in the Icelandic population corresponds to an absence of genetic divergence among individuals with different isotopic values and observed movement patterns.
The aim of our study was to investigate patterns of genetic structure among Icelandic killer whales. Tissue samples of photo-identified individuals allowed for the correlation of genetic data with existing data of photographic mark-recapture and individual isotopic niche.
Population subdivision based on microsatellite markers was estimated. Genetic diversity, mitochondrial DNA (mtDNA) haplotype frequency and individual data on observed movement patterns, social associations, and isotopic signatures were used to measure the divergence of putative genetic units. We interpret the identified genetic patterns among Icelandic killer whales, discuss potential ecological and behavioral processes driving them and how they correspond to the patterns of the best-studied Northeast Pacific ecotypes.

| Sample and data collection
Biopsy samples of wild Icelandic killer whales (N = 60) were collected from a research vessel in both herring overwintering-(winter) and spawning-(summer) grounds where killer whales can be seasonally found, coinciding with seasonal herring migration into those locations. In the winter, samples were collected off Grundarfjörður and blubber samples were collected from photo-identified individuals using a Remington rolling block system rifle ("Larsen" long-range biopsy system) with 35 or 40 mm sterilized biopsy tips in 2013 and an ARTS pneumatic darting system (Kvadsheim, Lam, Miller, Alves, & Others, 2009)

| DNA extraction and genetic sex identification
Total genomic DNA (N = 61) was extracted from skin samples using standard proteinase K digestion and phenol/chloroform methods (Sambrook, Fritsch, & Maniatis, 1989) modified for small samples by Baker et al. (1994). DNA was quantified on a NanoDrop ND-1000 Spectrophotometer and standardized to 10-20 ng/µl. The sex of individual whales was genetically identified using the protocol of Jayasankar, Anoop, and Rajagopalan (2008) which amplifies a 210-224 base pair (bp) fragment of the Y-chromosomespecific region (SRY) in males and a 442-445 bp size fragment of the ZFX/ZFY region in both sexes. Polymerase chain reaction (PCR) products were run on agarose gel stained with ethidium bromide (EtBr) and visualized under UV light.
For the smaller fragment, the quantities were the same except for the primers, which were increased to 0.3 μM of each. PCR amplifications were performed using a G-Storm GS1 thermal cycler (Gene Technologies) with an initial denaturation step at 95°C for 2 min, followed by a specific number of cycles (33 cycles for the longer fragments and 38 cycles for the smaller fragments) of denaturation at 95°C for 15 s, annealing for either 30 s at 54°C for the longer fragment or for 1 min at 57°C for the smaller fragment, and extension at 72°C for 1 min, followed by a final extension at 72°C for 5 min. Successful amplification was confirmed using agarose gel, EtBr staining, and UV visualization. Negative controls without DNA were included in all PCR plates to monitor for contamination during the PCR set up. Excess dNTPs and unincorporated primers from the completed amplifications were

| Microsatellite genotyping
Twenty-two loci were selected from the literature (Appendix S1).
Fifteen of these loci were previously used by Foote et al. (2011) for genetic differentiation of North Atlantic killer whales and seven by Parsons et al. (2013) for genetic differentiation among northern North Pacific killer whales. The microsatellites were arranged in five groups for multiplex PCR, according to the expected size range of each marker, the dye color used, and the annealing temperature in optimum PCR conditions (Appendix S1). We used QIAGEN® Multiplex PCR kits to amplify the loci with the fluorescent M13 tail single-reaction nested PCR method (Schuelke, 2000) with four different color-specific tails (Tysklind, 2009  for 90 s at 60°C for groups I-IV and 55°C for group V (Appendix S1), and extension at 72°C for 60 s, followed by 31 cycles of denaturation at 94°C for 30 s, annealing at 50°C for 90 s and extension at 72°C for 60 s, followed by a final extension at 60°C for 30 min. Successful amplification was confirmed by agarose gel, EtBr staining, and UV visualization. All amplifications included a negative control to detect contamination. Fragment analysis was conducted on an ABI 3730 DNA sequencer (Applied Biosystems) at Edinburgh Genomics (University of Edinburgh) using LIZ 500 (Applied Biosystems) as internal standard. Seven samples (>10% of the dataset) selected at random were re-amplified and re-genotyped for all loci to ensure consistent allele sizing. Alleles were sized using Peak Scanner Software 2 (Applied Biosystems) by two of the authors independently.
Allele frequencies were calculated in GENALEX v6.503 (Peakall & Smouse, 2006). All loci were inspected for scoring errors or null alleles using MICRO-CHECKER v2.2.3 (Van Oosterhout, Hutchinson, Wills, & Shipley, 2004). Deviations from Hardy-Weinberg equilibrium (HWE) and linkage equilibrium were tested using 1,000 iterations in GENEPOP v4.2 (Rousset, 2008). The analyses were conducted for the whole dataset and for each genetic unit identified by the clustering method (see Results section). Significance levels were corrected for multiple testing using the sequential Bonferroni correction (Holm, 1979) for these tests and all subsequent multiple comparisons in the study.

| Genetic differentiation and population subdivision
To identify genetic structure in the microsatellite markers' dataset, two different clustering methods were used. First, we used Bayesian clustering analysis for detection of genetically differentiated clusters (K) performed in STRUCTURE v2.3.4 (Pritchard, Stephens, & Donnelly, 2000). Ten independent runs for K values set from 1 to 10 were performed using a burn-in period of 100,000 iterations followed by 1,000,000 Markov chain Monte Carlo steps. The admixture model with correlated allele frequencies, recommended when population structure is likely subtle (Falush, Stephens, & Pritchard, 2003), was chosen, and no a priori information on the origin of the sample was indicated. The mean log-likelihood (L(K)) of each K, calculated in STRUCTURE HARVESTER v.0.6.94 (Earl & VonHoldt, 2012), was used as choice criterion to select likely number of K (Pritchard et al., 2000).
Secondly, we used a discriminant analyses of principal components method (DAPC) (Jombart, Devillard, & Balloux, 2010) in the package adegenet (Jombart, 2008) DAPC is a multivariate clustering method where individuals are clustered by genetic similarity not assuming any population genetic model, and efficiently detects genetic hierarchical structure (Jombart et al., 2010). We performed the DAPC analysis following the recommendations of Jombart and Collins (2015). Briefly, the most likely number of clusters was first assessed using a K-means method setting the maximum number of clusters to 40 and retaining all principal components. The most likely number was defined by the lowest BIC (Bayesian Information Criterion) value and the elbow in the BIC curve. Then, the genetic data were transformed using Principal Component Analyses and a linear discriminant analysis was performed on the retained principal components (no more than 80% were retained to avoid over-fitting). Each individual was assigned to a genetic unit according to its maximum membership probability.
The robustness of the division was confirmed by rerunning DAPC after removing one individual from the pairs of individuals showing a relatedness coefficient superior or equal to 0.45 as in Rosel, Hansen, and Hohn (2009) and Louis et al. (2014); pairwise relatedness values were estimated within each putative genetic unit (as in Louis et al., 2014, to ensure allele frequencies are affected by the inferred genetic structure) in KINGROUP v2_101202 (Konovalov, Manning, & Henshaw, 2004) using Queller and Goodnight's (1989) relatedness estimator in its symmetric form (Goodnight & Queller, 1999).
To characterize the level of differentiation of putative genetic units, we calculated pairwise and overall F ST values for microsatellite loci in FSTAT v2.9.3.2 (Goudet, 2001) using the computational methods of Weir and Cockerham (1984). The level of significance for pairwise F ST values was assessed using 3,000 permutations and analyses were also performed with the dataset excluding closely Microsatellite loci diversity indices were calculated for each putative genetic unit (and also per locus): the mean sample size per locus (n), the mean number of alleles per locus (k) and the observed and expected heterozygosities (H o and H e , respectively) in CERVUS (Kalinowski et al., 2007); the number of private alleles (PA) in GENALEX (Peakall & Smouse, 2006); the mean allelic richness (AR, i.e., mean number of alleles per locus averaged over the number of loci and adjusted for sample size) and the inbreeding coefficient (F IS ) in FSTAT (Goudet, 2001). Mean AR and H o values were compared among genetic units using two-sample Wilcoxon tests in R (R Core Team, 2015).
The distribution of mtDNA haplotypes within each unit was assessed and the relationship between haplotype and genetic unit was tested using a 2x3 Fisher's exact test in R (R Core Team, 2015).

| Genetic units and movement patterns, isotopic niche width and social segregation
The variation in movement patterns within each putative genetic unit was investigated by assigning an observed movement pattern to each sampled individual based on Samarra, Tavares, et al. (2017). Individuals were considered "seasonal" or "year-round" if photographed in Iceland only seasonally (either in the winter or in the summer) or if sighted in both herring overwintering-and spawning-grounds, respectively. We For the majority of the sampled individuals (N = 56), stable isotope ratios analyzed from skin samples (specifically of nitrogen denoted as δ 15 N and carbon denoted as δ 13 C) were available from Samarra, Vighi, et al. (2017). Differences in the isotopic niche width between putative genetic units were assessed using the area of standard ellipses corrected for sample size (SEA C ) and statistically tested by comparing the probability distributions of Bayesian estimates of SEA C (2,000,000 iterations and 10,000 burn-in iterations) with the SIAR (Parnell & Jackson, 2013) and SIBER packages (Jackson, Inger, Parnell, & Bearhop, 2011) (Csardi & Nepusz, 2006) in R (R Core Team., 2015).

| Sex-biased dispersal, recent migration rates, and changes in effective size
Sex-bias in dispersal among putative genetic units was tested using the biased dispersal option in FSTAT (Goudet, 2001). Differences in sex-specific F ST and variance of assignment index (vAI) between males and females from different genetic units, were tested by generating null distributions with 10,000 permutations. Since only adults or sub-adults were sampled, the whole dataset was used in this test. The mean number of successfully reproducing migrants per generation among movement patterns (N m ) was estimated using GENEPOP (Rousset, 2008). Recent migration rates (i.e., past 1-3 generations) among putative genetic units were estimated using the Bayesian method implemented in BAYESASS v3.0.4 (Wilson & Rannala, 2003). As recommended by Rannala (2015), preliminary runs were performed in BAYESASS to assess convergence and mixing using Tracer v1.6.0 (Rambaut, Drummond, Xie, Baele, & Suchard, 2018), which showed that convergence was usually reached before 100,000,000 iterations. Therefore, we performed five different runs with 300,000,000 Markov chain Monte Carlo iterations, a burn-in of 10,000,000 iterations, and a sampling interval of 2000 iterations.
Consistency in the results and effective convergence was confirmed for all runs.
Recent change in effective size of the putative genetic units was tested using BOTTLENECK v1.2.02 (Piry, Luikart, & Cornuet, 1999) using two different mutation models of microsatellite evolution: the infinite allele model and the two-phase model, which allows multiple-step mutations. Parameters were set to 10,000 repetitions, with 70% single-step mutations in the two-phase model and a variance of 12 for multiple-step mutations. Significance of results for both mutation models was assessed using the Wilcoxon test, which is more robust when used with few polymorphic loci (Piry et al., 1999).

| Sample collection and DNA extraction
The 61 photo-identified Icelandic killer whales were successfully genotyped and assigned a sex based on genetic analysis (18 females and 43 males): 31 individuals (nine females and 22 males) were seen in both seasonal herring grounds and the remainder of the individuals were seen only seasonally (11 were only seen in the winter, six were seen only in the winter and traveled to Scotland in the spring and summer, and 13 were only seen in the summer).

| mtDNA control region haplotype identification
Sequences of the control region (~611 bp) were generated for all individuals and two different haplotypes were identified, haplotypes 33 (published by Hoelzel et al., 2002) and 34 (published by Foote et al., 2009). The two haplotypes varied in only one site which was positioned in the longer fragment amplified in the mtDNA control region, with no variation seen in the smaller fragment. Genetic differentiation based on mtDNA was not calculated due to this low variation among samples.

| Microsatellite analysis
Of the 22 microsatellite markers, 19 were polymorphic in the Icelandic killer whale population (Appendix S1). The mean proportion of loci typed was 0.98 (nine samples failed to amplify at one locus, one at six loci and one at seven loci). All 61 genotypes were unique

| Genetic differentiation and population subdivision
The probability support produced by STRUCTURE was highest for K = 1, indicating no subdivision among samples (see Appendix S3).
However, three genetic units were identified by DAPC (Figure 2

| Genetic units and movement patterns, isotopic niche width and social segregation
Genetic unit 1 included 14 individuals seen year-round in Iceland and seven only seen seasonally (six in the summer and one in the winter).

Genetic unit 2 included 16 individuals seen year-round in Iceland and
nine only seen seasonally (seven in the winter and two in the summer).  social clusters had members assigned to both genetic units 1 and 3, four had members assigned to units 2 and 3, and two social clusters had members assigned to all of the three genetic units.

| Sex-biased dispersal, recent migration rates and changes in effective size
There was no significant sex-biased dispersal across observed movement patterns detected in FSTAT, since neither sex-specific F ST nor vAI values were significantly different between males and females (Table 2). Individual F IS values of each genetic unit were negative (Table 1) F I G U R E 5 Distribution of the estimates of standard ellipse area (‰ 2 ) based on 2,000,000 iterations for killer whales belonging to the three genetic units inferred by DAPC. The black dots represent the mode, the red crosses represent the area of the standard ellipse corrected for sample size (SEAc) of the real data, and the shaded boxes represent the 50%, 75% and 95% credible intervals (from dark to light gray, respectively) F I G U R E 6 Social associations among killer whales from each of the three genetic units. Nodes represent individuals and edges (links) represent existing associations between individuals. Thicker edges correspond to higher half-weight index (HWI) values.
Network was plotted using Fruchterman-Reingold force-directed layout (Fruchterman & Reingold, 1991)  correction for multiple tests). Furthermore, the allele frequency distribution had a normal L-shape for both genetic units, which is not expected if a recent bottleneck has occurred (Luikart, Allendorf, Cornuet, & Sherwin, 1998). Genetic unit 1 showed no significant indications of recent change in effective size (infinite alleles model: p-value = 0.04; two-phase model: p-value = 0.41; normal L-shaped allele frequency distribution).

| Fine-scale genetic structure
Our study supports weak fine-scale genetic differentiation among Icelandic killer whales into three genetic units. Population subdivision was supported by DAPC but not by STRUCTURE. STRUCTURE uses a model based-method that assumes Hardy-Weinberg equilibrium, which includes the strict assumption that populations are panmictic, discrete and reproductively isolated (Palsbøll, Peery, & Bérubé, 2010;Pritchard et al., 2000). When the assumptions are not met, the STRUCTURE analysis loses power to detect genetic structure. Furthermore, STRUCTURE can fail to identify genetic structure when: (a) allele frequencies vary gradually across a region (e.g., when there is isolation-by-distance) (Pritchard, Wen, & Falush, 2009); (b) sample size of genetic groups is limited or uneven (Kalinowski, 2011;Puechmaille, 2016;Waples & Gaggiotti, 2006); (c) gene flow is relatively high (Waples & Gaggiotti, 2006); and (d) mutation is low and genetic differentiation is limited (Almojil, Cliff, & Spaet, 2018;Barr et al., 2008;Latch, Dharmarajan, Glaubitz, & Rhodes, 2006;Waples & Gaggiotti, 2006). Contrarily, in cases of weak genetic differentiation, principal component analysis has proven to be a sensitive tool in detecting fine-scale structure (Novembre et al., 2008;O'Connor et al., 2015); for example, using DAPC, Almojil et al. (2018) found biologically meaningful genetic clusters where STRUCTURE failed to find any (as in this study). Therefore, we believe that for our dataset, DAPC was a more sensitive approach to accurately investigate finescale genetic differences among samples.
Overall, genetic units differed in isotopic signatures, mtDNA haplotype frequencies and observed movement patterns of the majority of its members. Although individuals associated more within the genetic unit (a result that can be affected by sampling bias of social groups), genetic units were not socially segregated and it is unlikely that association patterns were the main drivers of the genetic discontinuity observed. Since different social groups, as defined by Tavares et al. (2017), were sampled in this study and members of the same social clusters were assigned to different genetic units, it also seems unlikely that the sampling of some social groups more than others created the genetic structure observed among the samples.
The strongest level of genetic differentiation was found between individuals from genetic unit 3 and individuals from genetic units 1 and 2. Genetic unit 3 included individuals with overall higher δ 15 N values suggesting a broad diet (as discussed by Samarra, Vighi, et al., 2017), while individuals from genetic units 1 and 2 had a significantly narrower isotopic niche width and overall lower δ 15 N values, consistent with a diet predominantly composed of herring (as discussed by Samarra, Vighi, et al., 2017). Genetic units 1 and 2 had a significant predominance of mtDNA haplotype 33 and the majority of its individuals were seen year-round in Iceland. Contrarily, nearly all individuals in genetic unit 3 had mtDNA haplotype 34 and were seen only seasonally in Iceland; some of them are known to travel to Scotland where they were seen feeding upon marine mammals (Samarra & Foote, 2015;Samarra, Vighi, et al., 2017). It is important to acknowledge that since the time lapse during which the stable isotopic signal of a given diet remains in skin is only a few weeks (Browning, Dold, I-Fan, & Worthy, 2014;Giménez, Ramírez, Almunia, Forero, & de Stephanis, 2016), individual isotopic signatures might not represent the entire diet patterns of individuals. For example, the three individuals from genetic unit 3 (one seen in both seasons and two seen only in the summer) with δ 15 N values <13.5‰, the threshold encompassing all putative herring-specialists (Samarra, Vighi, et al., 2017), could still have a broad diet but they might have been feeding mainly on herring before they were sampled. However, it is also possible that not all individuals in genetic unit 3 have a broad diet and that there was insufficient data to further differentiate genetic unit 3 into separate subunits.
There was further significant genetic structure separating genetic units 1 and 2. Between these units, we found a small difference in the distribution of the δ 13 C values, suggesting some variation in geographical foraging area (Bearhop, Adams, Waldron, Fuller, & Macleod, 2004), When populations experience a recent reduction of effective size, allelic diversity is reduced faster than heterozygosity (Piry et al., 1999). However, our analysis found no convincing evidence of recent reduction in effective size of the genetic units. Other causes of excess of heterozygotes can be over-dominant selection (i.e., heterozygotes advantage due to higher fitness than homozygotes), outbreeding with a genetically distinct population, the presence of closely related or inbred family groups, the underlying social organizational level below that of the population and the composition of samples drawn from the real population (Beebee & Rowe, 2008;Kalinowski et al., 2007;Milkman, 1975;Parreira & Chikhi, 2015). We did not detect inbreeding, and the negative F IS values supported genetic unit exogamy (i.e., outbreeding) while a high number of migrants (N m ) indicated high gene flow.  Barrett-Lennard, 2000;Hoelzel et al., 2007Hoelzel et al., , 2002Parsons et al., 2013). The F ST values in our study are of similar magnitude to those in several mammal populations, assessed at the level of social groups with philopatry and gene flow among groups (Storz, 1999). Low genetic divergence has been commonly found in studies of marine species (e.g., Almojil et al., 2018;Vignaud et al., 2014;Ward, Woodwark, & Skibinski, 1994).

| Low genetic divergence
There was also extremely low mtDNA haplotype variation in our dataset, consistent with low worldwide mtDNA diversity patterns in killer whales (e.g., Hoelzel et al., 2007Hoelzel et al., , 2002Morin et al., 2010;Moura, Janse van Rensburg, et al., 2014). preferential mating within subgroups; (c) isolation-by-distance/time/ adaptation; and/or (d) periodic changes in connectivity among demes due to decadal changes in herring stock migration routes and separation into different herring feeding, spawning and wintering grounds.

| Adaptive divergence and stability of genetic structure
Delphinid population structure is likely driven by a combination of socioecological processes (e.g., Amaral et al., 2012;Fontaine et al., 2007;Gaspari, Azzellino, Airoldi, & Hoelzel, 2007;Louis et al., 2014;Möller, 2012;Moura et al., 2015Moura et al., , 2013. Ecological variation and resource specialization are key factors driving divergence in killer whale populations (Foote et al., 2016(Foote et al., , 2011Hoelzel et al., 1998). From our results, it is possible to speculate that the first-level of genetic division, separating genetic unit 3 from genetic units 1 and 2, corresponds to two main types of units: (a) individuals that potentially have a generalist diet and; (b) herring-specialists. Furthermore, the second-level of differentiation, separating genetic units 1 and 2, could be related to site-fidelity to herring grounds and subsequent distribution of the individuals along the coast of Iceland. For example, the genetic structure of southern right whales (Eubalaena australis) is maintained since site-fidelity to feeding areas is transmitted along matrilineal lines (Carroll et al., 2015;Valenzuela et al., 2009) et al., 2007, 2002), the Icelandic genetic units shared mtDNA haplotypes. Third, while between the mammal-eating and resident ecotypes and between resident subpopulations there is strong social avoidance (e.g., Ford et al., 1998;Ford et al., 2000), killer whales in Iceland assigned to different genetic units were not socially isolated and engaged in associations at least occasionally . Our results suggest that if killer whales in Iceland are progressing toward ecological divergence and niche specialization, this process is still at a very early stage. However, further exploration of different ecological niches by the different genetic units within the Icelandic population could eventually lead to adaptive variation and the formation of different ecotypes if it took place with stronger geographic and/or genetic isolation.
Our study raises the question of whether the fine-scale genetic structure observed among Icelandic killer whales is stable long-term and driven by local adaptation (e.g., through adaptation of foraging strategies for resource exploitation) or a subtle temporal structure driven by the contemporary distribution of fragmented seasonal herring grounds. Geographic overlap in Icelandic herring grounds and the fluidity of the feeding aggregations  may create a unique social opportunity that promotes gene flow between individuals from different genetic units, at least in seasonal herring grounds.
Gene flow dependent on seasonal spatiotemporal overlap of non-dispersing individuals is a pattern of mating seen within the Northeast Pacific ecotypes (e.g., Hoelzel et al., 2007;Parsons et al., 2013;Pilot, Dahlheim, & Hoelzel, 2010) and other low-dispersal species (e.g., longfinned pilot whales Globicephala melas, Amos, Schlötterer, & Tautz, 1993; African elephant Loxodonta africana, Nyakaana & Arctander, 1999). In this scenario, gene flow is not possible between individuals that never spatially overlap, for example, individuals that visit herring grounds at different times during the season or visit different herring grounds. Furthermore, herring is a prey known to change migration routes and abundance (Jakobsson & Stefánsson, 1999;Óskarsson, Gudmundsdottir, & Sigurdsson, 2009). Future changes in herring distribution might influence the spatiotemporal overlap or the level of dispersal from Iceland of different individuals and, consequently, the opportunities for social interactions and genetic admixture.
Metapopulation dynamics can arise when patches of a fragmented habitat are momentarily occupied and are followed by local extinction when conditions become unsuitable and then reoccupied again sometime later (Beebee & Rowe, 2008;McQuinn, 1997).
Continuous variation in genetic differentiation due to metapopulation dynamics would agree with the observed evolutionary patterns of several niche diversifications along genealogical lines of North Atlantic killer whales across thousands of years (Foote et al., 2013).
Additionally, metapopulation dynamics and local adaptation are not mutually exclusive and can simultaneously influence the structure of a population. For example, the sampling over consecutive years at multiple within-river sites of Atlantic salmon (Salmo salar) in the Sainte-Marguerite river (Canada) showed that the population had significant substructuring in both space and time (Garant et al., 2000). It is possible that a similar scenario is driving the patterns of differentiation at a relatively fine-scale within the Icelandic population. In such a situation, changes in movement patterns and geographic overlap could lead to future changes in the level of genetic divergence among Icelandic killer whales. Further studies investigating genetic differentiation at different time and space scales are needed to fully understand genetic structure among killer whales in Icelandic waters.

ACK N OWLED G M ENTS
We thank everyone that participated in the various field seasons,

CO N FLI C T O F I NTE R E S T
None declared.