Coast to coast: High genomic connectivity in North American scoters

Abstract Dispersal shapes demographic processes and therefore is fundamental to understanding biological, ecological, and evolutionary processes acting within populations. However, assessing population connectivity in scoters (Melanitta sp.) is challenging as these species have large spatial distributions that span remote landscapes, have varying nesting distributions (disjunct vs. continuous), exhibit unknown levels of dispersal, and vary in the timing of the formation of pair bonds (winter vs. fall/spring migration) that may influence the distribution of genetic diversity. Here, we used double‐digest restriction‐associated DNA sequence (ddRAD) and microsatellite genotype data to assess population structure within the three North American species of scoter (black scoter, M. americana; white‐winged scoter, M. deglandi; surf scoter, M. perspicillata), and between their European congeners (common scoter, M. nigra; velvet scoter, M. fusca). We uncovered no or weak genomic structure (ddRAD Φ ST < 0.019; microsatellite F ST < 0.004) within North America but high levels of structure among European congeners (ddRAD Φ ST > 0.155, microsatellite F ST > 0.086). The pattern of limited genomic structure within North America is shared with other sea duck species and is often attributed to male‐biased dispersal. Further, migratory tendencies (east vs. west) of female surf and white‐winged scoters in central Canada are known to vary across years, providing additional opportunities for intracontinental dispersal and a mechanism for the maintenance of genomic connectivity across North America. In contrast, the black scoter had relatively elevated levels of divergence between Alaska and Atlantic sites and a second genetic cluster found in Alaska at ddRAD loci was concordant with its disjunct breeding distribution suggestive of a dispersal barrier (behavioral or physical). Although scoter populations appear to be connected through a dispersal network, a small percentage (<4%) of ddRAD loci had elevated divergence which may be useful in linking areas (nesting, molting, staging, and wintering) throughout the annual cycle.

in general can be found nesting throughout most of northern Europe eastward through Russia to Siberia. Populations on both continents winter mainly in coastal regions: Atlantic and Pacific oceans along with the Great Lakes in North America; and primarily coastal areas of Black and Caspian seas with Asian white-winged scoters wintering in the northern Pacific Ocean (BirdLife International, 2018a, 2018b, 2018c. Scoters exhibit delayed breeding, low reproductive output with annually variable reproductive success, and are long-lived, which may make them more sensitive to factors that influence adult survival, such as environmental change (e.g., human disturbances or climate change; Anderson et al., 2015, Bordage & Savard, 2011, Brown & Fredrickson, 1997 and overharvest (Allendorf & Hard, 2009;Stott & Olson, 1972). Scoters are in apparent decline across their ranges, which has facilitated a Red List classification as vulnerable for velvet scoter (M. fusca) and near-threatened for black scoter (M. americana;BirdLife International, 2016a, 2016bBirdLife International, 2018a, 2018b, 2018c; status is unknown for common scoter (M. nigra). The influence of environmental change, harvest pressure, and other mechanisms regulating population dynamics, however, is poorly known for all scoter species. Scoter population recovery following declines, regardless of cause, may be longer when compared to other waterfowl (e.g., dabbling ducks) due to life-history characteristics of this group, and therefore, estimating levels of genetic structure among scoter populations is needed (Sea Duck Joint Venture, 2015).
Data regarding genetic connectivity will increase our understanding of species biology and provide critical information for predicting how these species may respond to future environmental and other disturbances.
Here, we present the first assessment of genomic connectivity within North American scoters to provide much needed insight into  (Table 1). Specifically, the black scoter exhibits a discontinuous range which may limit gene flow and result in genetic discontinuities. Conversely, the surf scoter has a continuous distribution suggestive of fewer inter-regional barriers to dispersal. Differences in the timing of mate selection have been hypothesized to account, at least in part, for differential patterns of genetic structure in geese (Ely & Scribner, 1994;Ely, Wilson, & Talbot, 2017;Wilson, Ely, & Talbot, 2018). Within scoters, black and surf scoter (and ducks in general) form pair bonds on the wintering grounds, which typically comprise individuals from multiple nesting locales, thus providing an opportunity for increased gene flow. Conversely, white-winged scoters pair during spring migration or soon after arrival at the nesting area (Brown & Fredrickson, 1997) and therefore may exhibit genetic structure among migratory pathways. This evidence suggests that scoters possess ecological and behavioral characteristics that have been shown to facilitate population structure as well as promote genetic connectivity across breeding regions in other waterfowl species (Sonsthagen, Talbot, Scribner, & McCracken, 2011;Wilson et al., 2018;Wilson, Gust, Petersen, & Talbot, 2016).
For this study, first, we used genome-wide scans (double-digest restriction-site-associated DNA sequences; ddRAD) and microsatellite genotype data to assess population genomic structure of the three North American scoter species across four regions (Alaska, Pacific, Central, and Atlantic). Second, we assessed levels of intercontinental dispersal between North American and European forms as vagrancy from migratory routes has been observed in all scoters (e.g., black scoter and surf scoter are observed in Europe, Scandinavia, and Russia; Anderson et al., 2015, Bordage & Savard, 2011. Third, as waterfowl are known for high hybridization rates (Gillham & Gillham, 1998;Kraus et al., 2012), we searched for evidence of potential hybridization or introgression among scoters.
Lastly, as information regarding annual linkages among areas is needed, we sought to identify high-resolution (i.e., divergent) loci with a strong potential for nesting area identification (e.g., Ruegg et al., 2014).

| Taxonomy and sampling
Of the five species of scoter (depending on taxonomic authorities), three occur in North America (black scoter, white-winged scoter, surf scoter) and two in Eurasia (common scoter, velvet/white-winged scoter). Traditionally, black scoter and common scoter were treated as conspecifics and were only recently regarded as different species, as they show distinct phenotypic characters (Chesser et al., 2016).
Similarly, taxonomic relationships between white-winged scoter and velvet scoter have varied over time; currently, they are considered conspecific by the American Ornithological Society and separate species by British Ornithological Union.
Samples were opportunistically collected as part of other research efforts, across the ranges of black, white-winged, and surf scoters in North America (Figure 2), as were representatives of the two European forms, common scoter (ddRAD N = 5;msats N = 19) and velvet/white-winged scoter (ddRAD N = 4; msats N = 20).
Samples were grouped into four broad North American regions: F I G U R E 2 Range (breeding: green and blue; wintering: orange) with location and number of samples assayed of black scoter (a, d), surf scoter (b, e), and whitewinged scoter (c, f)

| Library preparation and de-multiplexing
Genomic DNA was extracted from blood or tissue using a DNeasy Blood and Tissue kit following the manufacturer's protocols (Qiagen). Extractions for the ddRAD protocol were quantified using a NanoDrop 2000 Spectrophotometer (Thermo Fisher Scientific Inc.) to ensure a minimum concentration of 0.02 µg/µl and visualized on an 1% agarose gel for high molecular weight bands. Library preparation for multiplexing followed steps outlined in (DaCosta & Sorenson, 2014; also see Lavretsky et al., 2015Lavretsky et al., , 2016 Raw Illumina reads were de-multiplexed and processed using the computational pipeline described by DaCosta and Sorenson (2014; Python scripts available at http://github.com/BU-RAD-seq/ddRADseq-Pipeline) and following steps outlined in Lavretsky et al. (2015).
The pipeline clusters de-multiplexed and filtered reads into putative loci based on sequence similarity and genomic position as determined by BLAST, aligns reads within each putative locus, and infers genotyping for individual samples at each locus. Briefly, low-quality reads were filtered and identical reads collapsed for each sample.
Next, sequences were clustered into putative loci using the UCLUST function in USEARCH v. 5 (Edgar, 2010) with an -id setting of 0.85.
Clusters with identical or nearly identical BLAST hits (i.e., aligned to ± 50 bp on the same reference genome) were combined, which has been shown to minimize error associated with imposing a more TA B L E 1 Selected life-history attributes of North American scoters and hypothesized potential for population structuring. Species attributes that are unknown are denoted with a question mark  Harvey et al., 2015). Reads within each cluster (i.e., putative locus) were aligned using MUSCLE v. 3 (Edgar, 2004), and samples were genotyped using the Python script This information is produced in the "clustersummary.out" output file of the genotyping step of the DaCosta and Sorenson (2014) pipeline.
We used Geneious (Biomatters Inc.) to manually check and edit loci.
Doing so, allowed for the retention of many loci with insertions/deletions or high levels of polymorphism. Final datasets consisted of markers that had on average <10% missing genotypes.
Final output files (e.g., FASTA, NEXUS, and ADMIXTURE) were generated with custom python scripts that set a higher minimum sequencing depth to score an allele (Lavretsky et al., 2016). Specifically, to limit any biases due to sequencing error and/or allelic dropout, alleles with <5x coverage were scored as missing, such that a minimum of 10 reads was required to score a locus as heterozygous. Due to female heterogamy at sex chromosomes (Females = ZW, Males = ZZ), and our overarching goal to understand population structure among the species, all analyses were restricted to autosomal markers only.
Given that a higher number of loci are expected to be retained for within species than across species datasets, we aimed to test whether an increase in ddRAD markers provided additional resolution. Furthermore, because waterfowl readily hybridize across genera (Ottenburghs, Ydenberg, Hooft, Wieren, & Prins, 2015), combining all five species facilitates our ability to determine whether individuals of hybrid ancestry exist in our dataset. Consequently, a total of four ddRAD datasets were analyzed that included per-species alignments, as well as an alignment of all five species.

| Microsatellite genotyping
Genomic DNA was extracted following Medrano, Aasen, and Sharrow

| Population structure and diversity
Within and among species, we calculated composite pairwise esti- Pacific, Central, and Atlantic ( Figure 2). Finally, nucleotide diversity (π) was estimated in the R package PopGenome (Pfeifer et al., 2014) for chromosomally concatenated ddRAD autosomal loci.
Tests for HWE, LD, and F ST based on microsatellite data were corrected for multiple comparisons using Bonferroni correction (α = 0.05).
For ddRAD data, maximum-likelihood estimates of population as- were analyzed using an admixture model assuming correlated frequencies and sample location information as a prior with a 100,000 burn-in period, 1,000,000 Markov chain Monte Carlo iterations, and number of possible populations (K) ranging from 1 to 8; the analyses were repeated 10 times to ensure consistency across runs.
We followed the method of Evanno, Regnaut, and Goudet (2005) to determine the most likely number of clusters given the data.

| ddRAD dataset
For within-species analyses based on traditional species taxonomies (i.e., both North American and Eurasian forms of black scoter and white-winged scoter, and surf scoter), we recovered 3,487-3,999 ddRAD autosomal markers, comprising ~500K bp and ~15,000 SNPs (

| Population structure and diversity estimates
Within North America, pairwise estimates of genomic structure based on ddRAD loci were low across all species, but comparisons within the black scoter were higher relative to the other scoter species  F I G U R E 3 Average assignment probabilities of black/common scoter (a), surf scoter (b), and white-winged/velvet scoter (c) individuals from sampled regions (AK-Alaska; Pac-Pacific; Cen-Central; Atl-Atlantic; Eur-Europe) into two or three clusters inferred from ddRAD data in ADMIXTURE using all biallelic single nucleotide polymorphic sites along with pairwise estimates of genetic structure (ddRAD Φ ST below the diagonal; microsatellite F ST above the diagonal) among regions. Significant comparisons after Bonferroni correction (microsatellite data) are in bold text. Pairwise comparisons for white-winged scoter based on microsatellite data are not presented for nonbreeding and Alberta locales in the Pacific region  Table 2 for total biallelic SNPs per analysis). Analyses based on all biallelic SNPs and a random single SNP per locus yielded similar results    (Figure 3). In addition, black/common scoters and velvet/white-winged scoter individuals were assigned nearly exclusively to species-specific clusters ( Figure A3), although common scoter individuals clustered with black scoter individuals when all species were analyzed together ( Figure A4).
Calculated nucleotide diversity (Table 3) for black (π = 0.0031-0.0038) and white-winged (π = 0.0031-0.0036) scoters was slightly higher than those of surf scoters (π = 0.0028-0.0029). Diversity metrics were similar between the North American species and their European counterparts. Diversity metrics estimated from the microsatellite loci were similar across regions and species (Table 3).

| Among-species comparisons
The traditional species are highly differentiated from each other

| ddRAD dataset comparisons
Estimates of differentiation (R 2 = 0.99; p < 0.0001) and nucleotide diversity (R 2 = 0.98; p < 0.0001) were significantly correlated whether analyzing species separately (more loci) or together (fewer loci) ( Figure 3; Table 2; Figure A1). Moreover, when analyzing all samples together, ADMIXTURE assignments at a K of eight ( Figure A2) were concordant with those estimated in single species analyses (Figure 3; Figure A1). The latter results suggest that ADMIXTURE analyses are robust to differences in datasets.

| Comparisons within North America
The three North American scoter species exhibit weak or no discernible (ddRAD Φ ST < 0.019; microsatellite F ST < 0.004) genetic structure across their ranges, suggestive of high gene flow and connectivity among regions. The pattern of high connectivity among regions within North America is similar to patterns found in other sea duck species, such that spatial patterns in genetic variation are often not detected or have a weak signal at autosomal loci (Talbot, Sonsthagen, Pearce, & Scribner, 2015 Pacific wintering areas among years (Swoboda, 2007). A similar pattern is observed among king eiders that nest at Karrak Adventure lakes, Nunavut, Canada; the population is comprised of individuals that winter on the Atlantic and Pacific coasts, females switched (N = 6/20) winter areas between years (Mehl, Alisauskas, Hobson, & Kellett, 2004), and genetic structure was not detected among eastern and western populations (Pearce et al., 2004). Moreover, factors influencing an individual's migratory tendencies vary by species The lack of phylogeographic structure within white-winged scoters is particularly interesting because this species pairs during spring migration or summer, and therefore, we would expect some evidence of structure among regions as evident in other waterfowl species that share this characteristic (e.g., greater white-fronted goose, Wilson et al., 2018). Summer pairing would also be conducive of multiyear pairing and has been proposed (although not confirmed) for white-winged and surf scoters (Anderson et al., 2015;Brown & Fredrickson, 1997). Pair formation during spring migration or summer would restrict the availability of mates. However, recent studies using satellite telemetry found that two females used different spring migration routes between years (Meattey et al., 2018) and that male white-winged scoters migrated to different breeding areas between years (Sea Duck Joint Venture, 2015). If these behaviors are common in white-winged scoters, they may enable the formation of pair bonds between individuals from different breeding areas, serving to homogenize allele frequencies across nesting locales.
Population models generated to assess population declines in white-winged scoters nesting on Redberry Lake led the researchers to hypothesize that the population was rescued by female immigration (Alisauskas, Traylor, Swoboda, & Kehoe, 2004). A rescue (source-sink) dynamic via female dispersal mitigating population decline was also postulated to occur among nesting areas of spectacled eiders within the Yukon-Kuskokwim Delta, Alaska (Flint, Grand, Petersen, & Rockwell, 2016), and Chaun Delta, Chukotka (Solovyeva, Vartanyan, Frederiksen, & Fox, 2018). These studies indicate that females disperse among nesting areas (at least among areas that are experiencing decline), despite the presumption of high female philopatry within sea duck species (Eadie & Savard, 2015). The lack of genomic structure across North America indicates that intracontinental gene flow is occurring for both surf scoter and white-winged scoter. We assayed only putative autosomal loci and therefore cannot differentiate between male-or female-biased dispersal; data from maternally inherited markers are needed to confirm this hypothesis of female dispersal. Given that sea duck species share similar life-history characteristics, it seems likely that male-biased gene flow is also influencing the levels of genomic structure observed among regions within scoters.
Phylogeographic structure observed within black scoters at ddRAD loci is concordant with the disjunct breeding distribution. The barrier, whether behavioral or physical, between Alaska and Central/ Atlantic regions for black scoters appears to limit dispersal among areas as evidenced by the 17-fold higher Φ ST between segregated regions. Arctic and sub-Arctic taxa often exhibit a phylogeographic TA B L E 3 Indices of genetic diversity including the mean number of alleles (A), allelic richness (AR), and observed and expected heterozygosity (H o /H e ) based on 11 microsatellite loci, as well as nucleotide diversity (π) calculated using species-specific and concatenated ddRAD datasets (see Table 2) for black scoter, common scoter a , surf scoter, white-winged scoter, and velvet scoter individuals from sampled regions break located near the Mackenzie River, western Canada, which represents the eastern boundary of Beringia (Hewitt, 2004b); this signature is present, although subtle (ddRAD Φ ST = 0.017-0.019) within black scoter and not detected at microsatellite loci (Figure 3a).
Among sea duck species for which autosomal data are available (Pearce et al., 2004(Pearce et al., , 2005(Pearce et al., , 2008Peters et al., 2012;Scribner et al., 2001;Sonsthagen et al., 2011;Wilson et al., 2016), only the common eider exhibits a similar break in genetic structure among birds that winter in the Pacific versus Atlantic oceans (microsatellite F ST = 0.094, nuclear introns Φ ST = 0.098-0.120; Sonsthagen et al., 2011). Maintenance of genetic structure among regions in common eider was attributed to high female philopatry coupled with winter site fidelity, as partitions in the nuclear genome correspond to subspecific designations for the species .
While we only assayed autosomal loci and therefore cannot specifically hypothesize about female scoter dispersal tendencies, the spatial pattern of genomic diversity in black scoter suggests that the species may exhibit higher levels of winter site fidelity, which would decrease opportunities for pair-bonding with individuals from other nesting locations (i.e., lower incidence of intracontinental dispersal) than the other scoter species.
Finally, Holocene response following Pleistocene glacial cycling may also play a role in the shallow differentiation of populations of scoters. Most of northern North America was covered by the Cordilleran and Laurentide ice sheets through the last glacial maximum, and only recently colonized as species' distributions expanded into habitat made newly available by glacial retreat (Hewitt, 2004a). As scoters have only likely recently occupied northern North America (<11,000 years), insufficient time may have passed to accumulate regional level differences at nuclear markers (i.e., incomplete lineage sorting). Indeed, only a small percentage of loci (<3.9%) recovered exhibited elevated levels of divergence, which can be attributed to incomplete lineage sorting or dispersal. Both processes (dispersal and recent divergence), therefore, are likely playing a role in patterns of genetic diversity observed within scoters sampled across North America.

| Comparisons between continents
North American species are highly differentiated from their sons among locales with overlapping winter ranges were low and attributed to male-biased dispersal see also Scribner et al., 2001  . Genomic partitions are likely maintained by nonoverlapping (or nearly so) winter distributions (Collinson, Parkin, Knox, Sangster, & Helbig, 2006), which could limit opportunities to form pair bonds among North American and European varieties and ultimately restrict intercontinental dispersal. Examination of range-wide genomic structure among the other sea duck species (king eider, common goldeneye, long-tailed duck, and red-breasted merganser) is needed to confirm whether this pattern of high genomic structure between North American and Eurasian forms is a general characteristic of sea ducks, or unique to the scoters, common eider, and common merganser.

| Hybridization
Waterfowl are known for their propensity to hybridize with prezygotic mechanisms (e.g., courtship displays and vocalizations) maintaining species boundaries. Due to the high rate of genomic connectivity found within sympatric dabbling ducks (genus Anas), Kraus et al. (2012) coined the term "suprapopulation" where a group of species form a superspecies complex where natural hybridization occurs but without eroding species barriers. Scoter species are known to hybridize with each other; however, we did not detect any evidence of hybridization or introgression within or among traditional species (Figure 3), despite the traditional species pairs (whitewinged/velvet scoter and black/common scoter) demonstrating similar courtship displays (Collinson et al., 2006). Although our sample sizes for ddRAD analyses are small for the European forms, limiting our inferences, we failed to detect introgression among the three North American species based on our larger microsatellite genotype dataset (N = 457; see also Talbot et al., 2015). Courtship and copulation displays among black/common, velvet/white-winged, and surf scoters are diagnostic, perhaps providing a behavioral barrier to hybridization among the three traditionally accepted species (Collinson et al., 2006). Thus, although hybridization has been detected (white-winged scoter x surf scoter), our data support the supposition that these events are rare (Brown & Fredrickson, 1997) and do not appear to result in excessive introgression.

| Conservation implications
Scoters are migratory species that nest in remote regions of the Arctic, have experienced population declines, and are game species.
Furthermore, certain characteristics of scoter breeding biology (low reproductive output, delayed breeding, etc.) may limit their capacity to recover from population declines and events (either stochastic or deterministic) that reduce adult survival (Koneff et al., 2017).
Understanding linkages among key stages in the annual cycle are important to inform management strategies and conserve species, as events (i.e., disease, habitat quality, nutrition, and weather) during the nonbreeding season affect an individual's body condition, survival, and fecundity (Sedinger & Alisauskas, 2014). However, studies that investigate the level of population connectivity via dispersal in scoters (and waterfowl in general) must also contend with the lack of complete understanding about the relationship between migratory and dispersal behavior. Added to this layer of complexity is the fact that most waterfowl undertake a postbreeding remigial molt, where males, some nonbreeders, and failed female breeders migrate elsewhere to molt. Despite the lack of (or weak) genomic structure observed within North America that suggests that scoter populations are connected through a dispersal network that facilitates and/ or maintains panmixia, high-resolution loci (Φ ST > 0.20) among regions were uncovered (Figure 4). Analysis of these high-resolution loci assayed from nesting, molting, and wintering areas may provide an opportunity for researchers to further elucidate linkages among areas within the annual cycle (e.g., Ruegg et al., 2014) and provide insights on the composition of nesting areas among hunter-harvested individuals, enabling managers to develop harvest prescriptions that serve to reduce pressure on populations experiencing declines.

ACK N OWLED G M ENTS
We thank all of the researchers who provided samples through  Government.

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

AUTH O R S' CO NTR I B UTI O N
All authors conceived of the project, collected the data, performed the analyses, and contributed to the writing of the manuscript.

DATA ACCE SS I B I LIT Y
The data used in the present study are accessioned in NCBI F I G U R E A 1 Average assignment probabilities of black/common scoter, surf scoter, and white-winged/velvet scoter individuals from sampled regions (AK-Alaska; Pac-Pacific; Cen-Central; Atl-Atlantic; Eur-Europe) into one to three clusters inferred from 2,224 ddRAD loci analyzing a single randomly selected biallelic single nucleotide polymorphic site per locus in ADMIXTURE Alexander et al., 2009) F I G U R E A 2 Average assignment probabilities of black/common scoter, surf scoter, and white-winged/velvet scoter individuals from sampled regions (AK-Alaska; Pac-Pacific; Cen-Central; Atl-Atlantic; Eur-Europe) into three to eight clusters inferred from 2,224 ddRAD loci analyzing (a) all biallelic single nucleotide polymorphic sites (SNP) and (b) a single randomly selected biallelic single nucleotide polymorphic site per locus in ADMIXTURE Alexander et al., 2009  F I G U R E A 3 Average membership coefficient and ΔK (K = 2) of (a) black/ common scoter, (b) surf scoter, and (c) white-winged/velvet scoter individuals from sampled regions (Alaska, Pacific, Central, and Europe) inferred from 11 microsatellite loci in STRUCTURE (Pritchard et al., 2000). Individuals from the same sample locations are denoted with a black line Membership Coe cient F I G U R E A 4 Average membership coefficient, ΔK, and average likelihood given the data (LnP) of black/common scoter, surf scoter, and white-winged/ velvet scoter individuals from sampled regions (AK-Alaska; Pac-Pacific; Cen-Central; Atl-Atlantic; Eur-Europe) into three to eight clusters inferred from 11 microsatellite loci in STRUCTURE (Pritchard et al., 2000). The most likely number of clusters based on ΔK and LnP are denoted with an asterisk