Divergent selection and drift shape the genomes of two avian sister species spanning a saline–freshwater ecotone

Abstract The role of species divergence due to ecologically based divergent selection—or ecological speciation—in generating and maintaining biodiversity is a central question in evolutionary biology. Comparison of the genomes of phylogenetically related taxa spanning a selective habitat gradient enables discovery of divergent signatures of selection and thereby provides valuable insight into the role of divergent ecological selection in speciation. Tidal marsh ecosystems provide tractable opportunities for studying organisms' adaptations to selective pressures that underlie ecological divergence. Sharp environmental gradients across the saline–freshwater ecotone within tidal marshes present extreme adaptive challenges to terrestrial vertebrates. Here, we sequence 20 whole genomes of two avian sister species endemic to tidal marshes—the saltmarsh sparrow (Ammospiza caudacutus) and Nelson's sparrow (A. nelsoni)—to evaluate the influence of selective and demographic processes in shaping genome‐wide patterns of divergence. Genome‐wide divergence between these two recently diverged sister species was notably high (genome‐wide F ST = 0.32). Against a background of high genome‐wide divergence, regions of elevated divergence were widespread throughout the genome, as opposed to focused within islands of differentiation. These patterns may be the result of genetic drift resulting from past tidal march colonization events in conjunction with divergent selection to different environments. We identified several candidate genes that exhibited elevated divergence between saltmarsh and Nelson's sparrows, including genes linked to osmotic regulation, circadian rhythm, and plumage melanism—all putative candidates linked to adaptation to tidal marsh environments. These findings provide new insights into the roles of divergent selection and genetic drift in generating and maintaining biodiversity.


| INTRODUC TI ON
The role of differential adaptation to divergent selective environments in generating and maintaining biodiversity has become an increasing focus for evolutionary biologists over the past decade and has been termed ecological speciation (Funk, Nosil, & Etges, 2006;Rundle & Nosil, 2005;Schluter, 2001). While the concept of ecologically based divergent selection is not new (Mayr, 1942;Rundle & Nosil, 2005), disentangling the contribution of ecological forces from nonecologically based evolutionary forces (i.e., reductions in population size and genetic drift) remains a challenge. Despite these challenges, increasing accessibility and improvement of current sequencing technologies have allowed for the application of wholegenome sequencing to questions in natural populations (Campagna et al., 2017;Ellegren, 2014;Toews et al., 2016). The genomic era holds promise for ecological speciation research, as it has allowed for the detection of divergent signatures of selection on a genomewide scale. In recent years, the application of genomics to nonmodel systems has provided new insight into mechanisms underlying lineage-specific adaptations in a range of taxa (Cai et al., 2013;Li, Li, Jia, Caicedo, & Olsen, 2017;Liu et al., 2015;Qiu et al., 2012;Zhan et al., 2013). It has also provided insight into the genomic architecture of adaptation and speciation (Cruickshank & Hahn, 2014;Larson et al., 2017;Strasburg et al., 2012).
Under a scenario of ecological speciation, divergent selection will result in differential performance of individuals inhabiting alternative ecological niches (Arnegard et al., 2014;Nosil, 2012).
Reproductive isolation resulting from such adaptive divergence may occur either in sympatry, parapatry, or allopatry (Langerhans, Gifford, & Joseph, 2007;Nosil, 2007Nosil, , 2012. Comparison of the genomes of two ecologically divergent taxa spanning a selective habitat gradient provides valuable insight into the role of niche divergence in driving natural selection and speciation within a system. Characterizing genomic differentiation in recently diverged taxa is critical for increasing understanding of ecological speciation, as early acting barriers to gene flow have a larger effect on driving reproductive isolation than late-acting barriers (Coyne & Orr, 2004). In addition, identifying elevated regions of differentiation-and potential genomic regions under selection-is easier when baseline divergence is low, as expected for phylogenetically closely related and recently diverged taxa (Campagna et al., 2017;Ellegren et al., 2012;Poelstra et al., 2014;Toews et al., 2016).
Tidal marsh habitats in North America have undergone rapid changes since the last glacial maximum (Greenberg, 2006), and they are comprised of sharp ecotones that provide highly tractable opportunities for understanding underlying spatial patterns of genetic variation and adaptation. Because most tidal marsh endemics have colonized these habitats only after the rapid expansion of coastal marshes approximately 5,000-7,000 years ago (Malamud-Roam, Malamud-Roam, Watson, Collins, & Ingram, 2006), these systems provide opportunities to investigate patterns of recent and contemporary evolution. Environmental gradients across the saline-freshwater ecotone present extreme adaptive challenges to terrestrial vertebrates (Bayard & Elphick, 2011;Greenberg, 2006). Divergent selection among populations spanning these salinity gradients can be apparent in both physiological (i.e., pathways involved in osmotic regulation) and morphological (i.e., plumage-saltmarsh melanism, bill, and body size variation; Grenier & Greenberg, 2005;Grinnell, 1913) traits. These strong ecological gradients in tidal marshes provide a model system for applying comparative genomic analysis to investigate the role of ecological divergence in shaping species diversity.
Here, we investigated patterns of genome-wide differentiation between two recently diverged marsh endemics, the saltmarsh (Ammospiza caudacutus) and Nelson's (A. nelsoni) sparrow (~600,000 years; Rising & Avise, 1993). Although long considered a single species (AOU, 1931 A. caudacutus caudacutus-from southern Maine to New Jersey; and A. caudacutus diversus-from southern New Jersey to Virginia (Greenlaw & Rising, 1994). The prevailing evolutionary hypothesis (Greenlaw, 1993) suggests a history of vicariance for the saltmarsh and Nelson's sparrow, whereby an ancestral population spanning a coastal to interior range was split by Pleistocene glaciation, resulting in an isolated interior population. Following differentiation, this interior population then spread eastward back toward the Atlantic coast after recession of the Wisconsin ice mass, making secondary contact with ancestral coastal populations and establishing the current ranges and ecotypes within this species complex (Greenlaw, 1993).
Recent analyses of genetic and morphological characters indicate the strongest differences appear to correspond to habitat type, clustering the five subspecies into three groups: (a) the two freshwater, interior subspecies of Nelson's sparrow; (b) the brackish, coastal subspecies of Nelson's sparrow; and (c) the two saltwater, coastal subspecies of saltmarsh sparrow (Walsh et al., 2017).
While both species inhabit tidal marshes in sympatry, variation in habitat affinity, morphology, and behavior suggest a role for divergent selection and adaptation in this system. Specifically, the saltmarsh sparrow is a narrow niche specialist that has been associated with salt marshes over a longer evolutionary time frame (possibly 600,000 years, Chan, Hill, Maldonado, & Fleischer, 2006) compared with Nelson's sparrow, which uses a broader range of habitats including brackish and freshwater marshes and hayfields (Greenlaw, 1993;Nocera, Fitzgerald, Hanson, & Milton, 2007;Shriver, Hodgman, & Hanson, 2011). Further, Nelson's sparrow has lower nesting success in tidal, coastal marshes compared with saltmarsh sparrow, suggesting habitat-linked adaptive differences (Maxwell, 2018;Shriver, Vickery, Hodgman, & Gibbs, 2007;Walsh, Rowe, Olsen, Shriver, & Kovach, 2015). Due to these differences in niche specificity, differentiation of these sister species may have been largely driven by divergent natural selection.
Alternatively, changes in population size during vicariant isolation and colonization events may have increased the role of genetic drift in driving interspecific divergence. While previous work in this system has identified patterns consistent with divergent selection across the saline-freshwater ecotone (Walsh et al., 2017(Walsh et al., , 2015Walsh, Shriver, Olsen, & Kovach, 2016), the genome-wide pattern of differentiation between saltmarsh and Nelson's sparrows remains unknown. An understanding of the genomic landscape of these taxa will reveal the influence of demographic processes and divergent selection in ecological speciation.
We sequenced whole genomes of saltmarsh and Nelson's sparrows to investigate the role of divergent selection across an ecological gradient in shaping genome-wide patterns of divergence. We were interested in identifying genomic regions exhibiting elevated divergence due to selection across the saline-freshwater ecotone.
We predicted elevated divergence between saltmarsh and Nelson's sparrows in gene regions linked to known tidal marsh adaptations.

| Sample collection
For the reference genome, we sampled a male saltmarsh sparrow from the Marine Nature Center in Oceanside, New York, in July 2016. Blood was collected from the brachial vein and stored in Puregene lysis buffer (Gentra Systems, Minneapolis, MN). For genome resequencing, we sampled 20 individuals, 10 Nelson's sparrows (2 females and 8 males) and 10 saltmarsh sparrows (8 females and 2 individuals of unknown sex), from marshes along the northeastern coastline of the United States ( Figure 1; Table S1) during the breeding seasons (June-August) of 2008-2014. Nelson's sparrows were sampled from three populations in Maine. Saltmarsh sparrows were sampled from ten populations in Massachusetts, Rhode Island, Connecticut, New York, Maine, and New Hampshire. From each bird, F I G U R E 1 Breeding distributions and sampling locations. The breeding distributions of Nelson's and saltmarsh sparrows are shown in purple and orange, respectively. The hybrid zone along the New England coastline is shown in yellow. The breeding locations of Nelson's sparrows that were resequenced in this study are shown by the purple dots, the resequenced saltmarsh sparrows are shown by the orange dots, and the breeding location of the reference genome individual is shown in dark brown. The inset figures show the plumage of the two species we collected blood samples from the brachial vein and stored them on Whatman filter cards.

| Reference genome-library construction and assembly
For reference genome sequencing, DNA was extracted from blood stored in lysis buffer using a PureGene extraction kit (Gentra). Three libraries were constructed using an Illumina Nextera library preparation kit: one paired-end (PE) library with a 180 bp insert size and two mate-pair (MP) libraries with 3 and 8 kbp insert sizes. Each library was sequenced on a single HiSeq 2500 lane PE x 100 cycles by the Weill Cornell Medical College Core Genomics facility.
Genome assembly was performed with ALLPATHS-LG version 44849 (Gnerre et al., 2011;Ribeiro et al., 2012) using the default parameters. ALLPATHS-LG takes raw data as input, without prior adapter removal and trimming. We used bioanalyzer results to estimate the insert size and expected standard deviation, which are required input for ALLPATHS-LG. The assembly was completed in six days on a 64-core computer (1,024 GB RAM, 19 TB hard drive) from the Cornell Computational Biology Service Unit BioHPC Lab.
We obtained assembly statistics with Quast version 2.3.

| Reference genome-annotation
We annotated the saltmarsh sparrow assembly with the MAKER pipeline v 2.31.9 (Cantarel et al., 2008). Gene models were cre- University. Transcriptome assembly was performed on the combined datasets with CLC Genomics Workbench (V5.1.2.) using CLC Assembly Cell 4.0 set to default parameters. Genes were predicted with SNAP v2013-11-29, using an iterative training process inside MAKER v2.31.9, and Augustus v3.2.2_4, using a hidden Markov model from the chicken. This produced a total of 15,414 gene models, which included 85.2% complete BUSCO v2.0 (Simão, Waterhouse, Ioannidis, Kriventseva, & Zdobnov, 2015) genes and a further 9.3% which were fragmented, when assessed against the aves_odb9 lineage data and specifying the white-throated sparrow (Zonotrichia albicollis) as the most closely related species in the database. The saltmarsh sparrow assembly was aligned to the zebra finch genome using Satsuma Chromosembler (Grabherr et al., 2011). Using this approach, scaffolds were mapped to chromosome coordinates via synteny. Thus, scaffolds were assigned to chromosomes for removal of sex-linked markers; however, population genetic analyses were conducted at the scaffold level.

| Whole-genome resequencing and variant discovery
We sequenced the genomes of an additional 20 individuals, 10 saltmarsh sparrows and 10 Nelson's sparrows. DNA was extracted ac.uk/proje cts/fastqc).
The average alignment rate was 97.3 ± 0.90 and 97.8 ± 0.43 across all saltmarsh and Nelson's sparrow individuals, respectively. After aligning sequences to the saltmarsh sparrow reference genome, depth of coverage ranged from 5 to 34X (18.7 and 34X for the two individuals sequenced in single lanes and 5-8X for the remainder that were multiplexed for sequencing; Table S1).
BAM files were sorted and indexed using Samtools version 1.3 , and PCR duplicates were marked with Picard Tools version 2.1.1 (http://picard.sourc eforge.net). We realigned around indels and fixed mate pairs using GATK version 3.5 (McKenna et al., 2010). SNP variant discovery and genotyping for the 20 resequenced individuals were performed using the unified genotyper module in GATK. We used the following filtering parameters to remove variants: QD < 2, FS > 40.0, MQ < 20.0, and HaplotypeScore > 12.0. In addition, variants that were not biallelic, had minor allele frequencies <5%, mean coverage less than 3X or with coverage greater than two standard deviations above the mean, and more than 20% missing data across all individuals were filtered from the data set. This resulted in 7,240,443 SNPs across the genome; the mean SNP depth averaged over individuals was 10.7 (standard deviation: 11.2). Using mapping results from Satsuma, we removed all SNPs located on the Z chromosome to avoid any bias that may be introduced by analyzing a mix of male and female individuals, resulting in a final dataset containing 6,256,980 SNPs.

| Population genomics
Principle component analysis (PCA) was performed on all SNPs using the SNPRelate package in R (R Development Core Team, 2017). We identified divergent regions of the genome by calculating F ST values for nonoverlapping 100 kb using VCFtools version 0.1.14 (Danecek et al., 2011). Descriptive statistics for the 100 kb windows, including the number of fixed SNPs, the proportion of fixed SNPs in coding regions, nucleotide diversity (π), Tajima's D, mean observed heterozygosity (H obs ), and mean expected heterozygosity (H exp ), were calculated using VCFtools and R. We calculated absolute divergence (D xy ) using a custom python script (S. Martin, https ://github.com/simon hmart in/genom ics_general). Divergent peaks were visualized using Manhattan plots, which were constructed using the R package qqman. We discarded regions with less than two windows and windows with less than 10 SNPs.

| Characterizing divergence between Nelson's and saltmarsh sparrows
To fully assess genome-wide patterns of differentiation between saltmarsh and Nelson's sparrows and to identify potential genes that play a role in adaptive divergence, we employed a multistep approach to identifying and characterizing candidate regions of interest. First, we defined candidate genomic regions under selection if they contained window-based F ST estimates above the 99th percentile of the empirical distribution. Second, we estimated the density of fixed SNPs (df) in the same 100-kb windows following Ellegren et al. (2012) and identified windows above the 99th percentile of the df distribution. For both approaches, elevated windows were inspected in Genious version 9.1.5 (Kearse et al., 2012). We compiled a list of gene models within 50 kb of each elevated region and obtained information on these annotations from the UniProt database (http://www.unipr ot.org/).
We performed GO analyses of divergent windows (Table 1) using the Web-based GOfinch tool (http://bioin forma tics.iah.ac.uk/tools/ Gofinch). Lastly, we compiled a list of candidate genes hypothesized to be important for tidal marsh adaptations, including genes linked to reproductive timing (circadian rhythm genes), osmotic regulation, salt marsh melanism, and bill morphology. Genes were chosen based on previous research done in this system (Walsh et al., 2017(Walsh et al., , 2016) and a literature review. We estimated F ST and df for each candidate gene of interest plus 20 kb upstream and downstream of the gene, and compared our divergence estimates to the genome-wide average.

| Modeling demographic history
An investigation of population size history was performed with the pairwise sequentially Markovian coalescent (PSMC) model (Li & Durbin, 2011). We used a single individual from each species, each with mean coverage ≥18X and missing data < 1%. Coverage ≥ 18X and <20% missing data have been shown to be important for accurate inference using PSMC (Nadachowska-Brzyska, Burri, Smeds, & Ellegren, 2016). Samtools and bcftools were used to create a diploid consensus genome for each individual. For the PSMC inference, we used time intervals denoted by −p "60*1 + 2 + 2" meaning that the first sixty atomic time intervals are of length 1, followed by two intervals with length 2. All other parameters were kept as defaults. We ensured that at least 10 coalescent events occurred in each interval after 20 iterations and performed 100 bootstrap replicates. For plotting, we used a mutation rate of 7.59 × 10 −9 subs/site/generation and a generation time of 2.3 years. The mutation rate was calculated using the mean mutation rate of passeriformes, which was estimated to be 3.3 × 10 −3 substitutions per site per million years (Zhang et al., 2014). Generation time was estimated using the equation 1/m + b where m is mortality, estimated to be approximately 0.5 for adults and one and a half times that for juveniles, and b, the age of first breeding, is equal to 1 (Greenlaw, Elphick, Post, & Rising, 2018).
To investigate whether the changes in effective population size seen in the PSMC were real or due to admixture between the two species, we used approximate Bayesian computation, implemented in DIYABC v2.1.0 (Cornuet et al., 2014) to estimate the posterior probabilities of three different scenarios for the demographic history of the populations. DIYABC has been shown to have limited ability to discriminate between a large number of complex models (Cabrera & Palsbøll, 2017) so we limited our investigations to three, relatively simple scenarios ( Figure 2). The scenarios that we tested were as follows: (a) a "simple split" scenario, similar to Greenlaw's (1993) hypothesis, where the ancestor to both species was subdivided into an inland, ancestral Nelson's population and a coastal, saltmarsh sparrow population, Nelson's sparrow population further splitting, giving rise to the inland A. n. nelsoni and the coastal, sampled A. n. subvirgatus; (b) a "change in Nelson's Ne" scenario which follows the simple split scenario except, after diverging from A. n. nelsoni, A. n. subvirgatus goes through an increase and then a decrease in effective population size as seen in the PSMC plot (see Results); and finally (c) an "admixture" scenario which also follows the simple split model but after A. n. subvirgatus has diverged from A. n. nelsoni, an admixture event between A. n. subvirgatus and the saltmarsh sparrow creates a new, admixed lineage of A. n. subvirgatus. In all scenarios, A. n. nelsoni is a ghost lineage ( Figure 2) since it is not sampled in this study.
Three sets of 2000 autosomal SNPs were selected, and the DIYABC simulations were repeated with each set to ensure results were robust to the subset of SNPs chosen. High-quality SNPs were selected that had no missing data, had mean depth >8 and <12, were sequenced to a depth of at least five in all individuals (high-confidence genotype calls are necessary for demographic history reconstruction), were within 1 standard deviation of the mean F ST to ensure they were effectively neutral, were at least 20 kb apart to reduce linkage, and then chosen at random without replacement.  Wiest et al., 2016, Erskine, 1992, Rivard, Shaffer, & Falardeau, 2006 and the results of test simulations, which indicated that large priors were needed for these parameters. We had little information about the effective population size of the ancestral Nelson's sparrow (N AncestralNESP ) so we used a wide prior of 1,000-150,000 individuals.
PSMC indicated that the ancestor to saltmarsh and Nelson's sparrows (N Ancestral ) may have had very large effective population size, so we used 100,000-500,000 as a prior range for this parameter.
In the "change in Nelson's N e " scenario, we used the PSMC plot

| RE SULTS
The final reference genome assembly generated by ALLPATHS-LG consisted of 44,080 contigs with an N50 of 66.4 kb and 2,672 scaffolds with an N50 of 8.427 Mb. Contig length was 1.03 Gb, and the total scaffold length was 1.07 Gb. Based on a 1.0 Gb genome, sequencing coverage for the assembly was 83X. Statistics for the final assembly are included in Table S2. We assessed the completeness of our reference assembly by searching for a vertebrate set of 3,023 single-copy orthologs using BUSCO version 1.2 (Simão et al., 2015).
Our final saltmarsh sparrow reference genome contained a single and complete copy of 83.6% of the genes in the vertebrate set and a partial copy of an additional 7.5% of the genes in this set. We found 0.5% of the BUSCO vertebrate genes more than once within the reference genome, and 8.8% of the BUSCO genes were missing from the saltmarsh sparrow reference.
Populations of saltmarsh and Nelson's sparrows exhibited high levels of genome-wide divergence. Individual sparrows clustered strongly by species in a PCA based on 7.2 million SNPs, with PC axis one explaining 40.9% of the observed variation ( Figure S1). Genomewide average F ST (across all autosomal SNPs) was 0.32 (Figure 3a;  (Tables S3 and S4). We also identified 94 windows across 42 scaffolds exhibiting values of df greater than the 99th percentile of the overall distribution (Tables S5   and S6). When combining these two criteria, we found 33 windows across 16 scaffolds that showed elevated divergence using both F ST estimates and the density of fixed SNPs (hereafter elevated regions, Table   S7). When comparing elevated regions to the rest of the genome, we  Table 2). Within these elevated regions, we identified 19 genes with putative roles in adaptive differences between the species (Table 1, Figure 4). Of these genes, nine were linked to osmoregulatory function, two were linked to reproductive differences between the species, two were linked to immune response, and one was linked to circadian rhythm. Enrichment analyses of these genes identified several pathways, including some linked to a priori hypotheses: intracellular calcium-activated chloride channel activity (potentially linked to osmoregulatory function; p = .0067), regulation of circadian rhythm (potentially linked to nest initiation in relation to tidal cycles; p = .023), and sodium:potassium-exchanging ATPase activity (potential link to osmoregulatory function; p = .023). Additional genes with potential tidal marsh or other adaptive functions were identified by either the F ST or df criteria alone (Tables S3-S6). Nucleotide diversity was decreased within these elevated windows, although not statistically significantly given the high genome-wide variation (Table 2).
This could be evidence for selective sweeps in these regions. Lastly, based on a literature review and a priori predictions, we identified several a priori candidate genes linked to putative tidal marsh adaptations (Table 3). Of these candidates, only two had an F ST value greater than the upper bound of the 95% confidence interval for genome-wide F ST (Table 3) Analyses of population size history with PSMC revealed that both species diverged from an ancestral population that had much larger effective size (150,000-200,000) than either extant species ( Figure 5).
The saltmarsh sparrow appears to have had a long history of relatively an admixture event, occurred at this time, followed by a decrease in population size to approximately 20,000 individuals 10,000 YA. Our DIYABC analysis confirmed that the hump in the PSMC plot was likely an admixture event. The "admixture scenario" was highly favored (posterior probability = 1.000) over the "change in Nelson's Ne" scenario (posterior probability = 0) and the "simple split" scenario (posterior probability = 0). Confidence in scenario choice was high (posterior predictive error of scenario choice = 0.000), and this was true for all three sets of 2000 SNPs that were used. However, the timing of the admixture event differed for the DIYABC analysis, which estimated the event to have occurred within the last 10,000 years (Table 4), rather than 50,000-100,000 YA as indicated in the PSMC plot.

The estimates of time since divergence of the saltmarsh and
Nelson's sparrow lineages were similar for the two analyses, with PSMC suggesting divergence at ~1 million YA and the three median estimates from DIYABC ranging 478,400-655,500 YA (95% HPD range of 377,200-1,281,100 YA; Table 4). Note that the lowest estimate stemming from the 3rd SNP subset may not be robust, given that the mean was approaching the lower bound of the prior. The DIYABC parameter estimates for the saltmarsh sparrow, Nelson's sparrow, and ancestral Ne were relatively consistent across the three SNP subsets (Table 4) and largely consistent with the estimates from PSMC, with Nelson's sparrow having a smaller population size than saltmarsh sparrow, although note that DIYABC estimates are averaged across the full time span of the analysis, rather than for a particular time period, as for PSMC. Median saltmarsh sparrow Ne averaged 110,000 (95% HPD range 34,900-149,000) and median Nelson's sparrow averaged 23,600 (95% HPD range 4320-48,400) across the three SNP subset analyses. Ancestral Ne estimates averaged 391,000 (95% HPD range 121,00-490,000). One of the SNP subsets for saltmarsh sparrow Ne and two of the SNP subsets for ancestral Ne produced means that were approaching the upper bound of the prior.

| D ISCUSS I ON
Whole-genome comparisons of saltmarsh and Nelson's sparrows revealed high baseline divergence between species (genome-wide F ST = 0.32). We identified a high density of fixed SNPs (~234,000),

TA B L E 3 (Continued)
F I G U R E 5 Demographic history inferred using PSMC. We performed 100 bootstrap replicates (thin lines) for each species. Saltmarsh sparrows (orange) appear to have had a long history of relatively small population size since they diverged from the common ancestor, which had much larger effective population size. Nelson's sparrows (purple) have a more complex history, with either a change in population size or a change in population structure occurring approximately 50-100 kya which appear to be uniformly distributed across the genome.
These patterns share commonalities with that observed in the collared (Ficedula albicollis) and pied (F. hypoleuca) flycatcher, which exhibit high baseline divergence in allopatry (F ST = 0.357 ;Ellegren et al., 2012), a large number of fixed SNPs (239,745), low rates of heterospecific pairings, and fitness reductions in hybrids (Svedin, Wiley, Veen, Gustafsson, & Qvarnström, 2008). The latter two observations have also been made in saltmarsh and Nelson's sparrows (Maxwell, 2018;Walsh, Maxwell, & Kovach, 2018;Walsh et al., 2016). Conversely, flycatchers may exhibit a deeper divergence time (<2 million years; Ellegren et al., 2012) than that estimated for saltmarsh and Nelson's sparrows-500,000-1 million years based on our DIYABC and PSMC analyses in this study, consistent with 600,000 years based on mitochondrial differentiation (Rising & Avise, 1993). Our findings suggest higher genome-wide patterns of divergence than between other shallowly diverged, hybridizing taxa. For instance, Poelstra et al. (2014) and Toews et al. (2016) found only a small number of elevated regions and a small num-  May not be well estimated from this SNP subset, as the mean was approaching the lower bound of the prior.
splitting, or founder events when re-colonizing tidal marsh habitats) may have provided a prominent role for drift in shaping genomic divergence between these taxa.
In addition to historical processes, patterns of contemporary gene flow may also shape the genomic landscape. Despite a 200-km hybrid zone between saltmarsh and Nelson's sparrows (Hodgman, Shriver, & Vickery, 2002;Shriver et al., 2005), our sampled individuals are predominantly from allopatric populations separated by both geographic distance and a habitat gradient that may pose a selective barrier in this system (Greenlaw, 1993;Walsh et al., 2015). Divergence via drift in allopatry therefore may have progressed unfettered from potential homogenizing effects of gene flow, resulting in the observed patterns of high divergence across the genomic landscape (sensu Feder, Egan, & Nosil, 2012). Support for this idea comes from a study comparing whole genomes of saltmarsh and Nelson's sparrows from allopatric and sympatric populations (Walsh, Kovach, Olsen, Shriver, & Lovette, 2018). That study found that only about 5% of the fixed differences found in allopatry are present in sympatric populations, suggesting that when populations co-occur contemporary gene flow homogenizes all but a small portion of the genomic landscape, which likely comprises important barrier loci (loci important in reproductive isolation ;Feder et al., 2012;. Thus, drift-related to both historical and contemporary processes-seems to be an important factor in shaping the landscape of genomic divergence we have observed in these Ammospiza sparrows. Despite a relatively high baseline of genome-wide divergence between saltmarsh and Nelson's sparrows, we detected putative candidate regions of adaptation housed within elevated regions of differentiation (expressed both as elevated estimates of F ST and as concentrated densities of fixed SNPs). This suggests that demographic processes alone are not responsible for the observed genomic landscape and supports a role for divergent selection in shaping the observed patterns of genome-wide divergence between these taxa. Agreement on whether ecological divergence should be manifested in localized versus genome-wide differentiation is generally lacking, and it is likely that different processes are operating at different stages of the speciation continuum Hemmer-Hansen et al., 2013;Via, 2012). Genomic islands of divergence-discrete regions of the genome with elevated divergence harboring clusters of loci underlying ecological adaptations-have been found in ecotypes or sister species (e.g., Hohenlohe, Bassham, Currey, & Cresko, 2012;Jones et al., 2012;Larson et al., 2017;Nosil, Harmon, & Seehausen, 2009). These islands are commonly thought to arise through divergence hitchhiking, whereby strong selection in conjunction with reduced recombination (due to linkage disequilibrium) results in coordinated evolution of multiple genes in the same genome region (Via, 2012), although other mechanisms, including regions of low genetic diversity, may also explain their presence (Cruickshank & Hahn, 2014). The prevalence, functional significance, and underlying mechanisms of these islands of divergence are not yet fully understood (Cruickshank & Hahn, 2014;Larson et al., 2017), and adaptive loci have also been found to be distributed more evenly across the genome (Strasburg et al., 2012). The pattern of high background differentiation coupled with peaks of elevated divergence distributed throughout the genomes of saltmarsh and Nelson's sparrows is consistent with expectations for populations that diverged in allopatry, in contrast to the clustered patterns of divergence found in species that diverged with ongoing gene flow (reviewed in Harrison & Larson, 2016). Under a scenario of ecological speciation, the differentiated genes within these elevated genomic regions underlie adaptations to the differential selective pressures faced by saltmarsh and Nelson's sparrows in allopatry.
We hypothesize that several of the observed genetic differences represent ecologically favored alleles in either saline or freshwater habitats, alleles underlying adaptive behaviors to tidal versus nontidal environments, or alleles that represent other ecological differences between the species.
Using a multistep approach, we identified strong candidate regions for tidal marsh adaptations between saltmarsh and Nelson's sparrows.
We identified several candidates with elevated F ST estimates or numbers of fixed differences, with our most compelling candidates exhib- been previously found to exhibit reduced introgression and increased selection across the saltmarsh-Nelson's sparrow hybrid zone (Walsh et al., 2016), offering further support for an adaptive role for this gene.
ALB is linked to the regulation of osmotic pressure of blood plasma and is a gene that was found to be under positive selection in the saline tolerant crab-eating frog (Fejervarya cancrivora) compared with a morphologically similar, saline-intolerant sister species (Shao et al., 2015).
Another candidate gene linked to tidal marsh adaptation identified in this study was CRY1, which is an important component of the circadian clock. In salt marshes, flooding affects nests during the highest spring tides; during this time, the entire marsh is flooded and nests can be inundated with water for an hour or two (Gjerdrum, Sullivan-Wiley, King, Rubega, & Elphick, 2008). Monthly tidal flooding, therefore, is the leading cause of nest failure and an important driver of overall population trajectories for tidal marsh sparrows (Greenlaw & Rising, 1994;Shriver et al., 2007); synchronizing the 23-to 26-day nesting cycle with the approximate 28-day monthly tidal cycle is critical for individual fitness. Saltmarsh sparrows have greater nesting synchrony with tidal cycles compared with Nelson's sparrows, and this synchrony is associated with increased nesting success (Shriver et al., 2007;Walsh et al., 2016). Biological clocks are important for ensuring that the passage of time is synchronized with periodic environmental events (Kumar et al., 2010); as such, it is reasonable to hypothesize that CRY1 may play a role in the divergent nesting synchrony observed between these species. Expression patterns of CRY1 were shown to fluctuate with lunar periodicity in a lunar-synchronized spawner, the golden-lined spinefoot (Siganus guttatus; Ikegami, Takeuchi, Hur, & Takemura, 2014), supporting a link between CRY1 and lunar cycles.
The final putative candidate for tidal marsh adaptation that we identified was TYRP1, which is an important component of the melanin biosynthesis pathway and has been found to play an important role in determining plumage color in birds (Bourgeois et al., 2016;Minvielle, Cecchi, Passamonti, Gourichon, & Renieri, 2009;Xu, Zhang, & Pang, 2013). This finding supports previous work identifying SLC45A2 (another gene associated with plumage melanism) to be under divergent selection in this system (Walsh et al., 2016).
Saltmarsh and Nelson's sparrows have subtle plumage differenceswith saltmarsh sparrows showing darker breast and flank streaking and face coloration than Nelson's sparrows (Shriver et al., 2007)consistent with phenotypic patterns observed in other vertebrates spanning tidal marsh gradients (Grinnell, 1913;Luttrell et al., 2014).
Tidal marsh taxa are grayer or blacker than their upland relatives, due to a greater expression of eumelanin relative to pheomelanin (Greenberg & Droege, 1990). This salt marsh melanism is thought to confer adaptive benefits to tidal marsh birds either via enhanced predator avoidance (from background matching with the gray-black marsh mud; Greenberg & Droege, 1990) or resistance to feather degradation (increased melanin concentrations slow degradation rates by salt-tolerant feather bacteria; Olsen, Greenberg, Liu, Felch, & Walters, 2009;Peele, Burtt, Schroeder, & Greenberg, 2009).
In addition to candidate genes that were in line with our a priori predictions for tidal marsh adaptations, we also identified several genes under selection that are related to spermatogenesis.
Differences in mating strategies in saltmarsh and Nelson's sparrows (Greenlaw, 1993) coupled with strong assortative mating in the hybrid zone (Walsh et al. in press) support a role for premating barriers, which could include sperm competition or incompatibilities between the species. These candidate regions, coupled with the osmoregulatory and circadian rhythm candidates described above, will provide directions for future research in this system, including candidate gene work in more individuals across habitat types. Further, the full suite of genes identified by either F ST or df criteria alone (Tables S3-S6) provides additional putative candidates for further investigation.
Our results demonstrate the effects of both ecological divergence and drift in driving high baseline levels of differentiation between two closely related sister taxa. The candidate genes for osmoregulatory response, circadian rhythm, and melanistic plumage that we identified in this study likely represent important lineage-specific adaptations and add to a body of empirical evidence describing underlying mechanisms driving adaptation to harsh environments (Tong, Tian, & Zhao, 2017;Wu et al., 2014;Yang et al., 2016). Our findings, in particular, contribute to a growing list of candidate genes linked to salt tolerance and osmoregulation (Ferchaud et al., 2014;Gibbons, Metzger, Healy, & Schulte, 2017;Kahle, Rinehart, & Lifton, 2010;Walsh et al., 2019), suggesting that modifications to several pathways can result in adaptations to saline environments.
Demographic analyses underscore the influence of population history, including multiple admixture events and population declines on the genomic landscape. These analyses also highlight that although the two species have persisted with relatively small effective population sizes since their initial divergence from a common ancestor, those population sizes were at least an order of magnitude higher than those of saltmarsh and Nelson's sparrows today (Wiest et al., 2016). High genome-wide divergence and a high proportion of fixed SNPs distributed across the genomes of saltmarsh and Nelson's sparrows offer new perspectives into processes shaping the genomic landscape and offer empirical evidence for the shared roles of ecological divergence and demography in shaping evolutionary processes and genetic variation between these taxa.

ACK N OWLED G M ENTS
We thank Kazufusa Okamoto for transcriptome assembly and Stephen Simpson for library development of resequenced genomes.

CO N FLI C T O F I NTE R E S T
We have no competing interests.

AUTH O R CO NTR I B UTI O N S
A.I.K. conceived and designed this study, with input from J.W.,