Discussion

We report on the first investigation of population genetic diversity and structure among all three spawning ecotypes of pike (freshwater, anadromous and brackish water resident) using information based on both neutral and adaptive molecular markers (RADseq SNPs). We document patterns of genetic structure at different spatial scales and identify the eco-evolutionary drivers of the genetic differentiation among and within ecotypes. Besides providing rare evidence of contrasting patterns of neutral and adaptive genetic structure, results exemplify how separate analyses of coding and non-coding variation can help disentangle the complex interplay of different stochastic and deterministic contributing processes. Specifically, we found that for neutral differentiation allopatry seems to play an important role. For the sympatric Baltic Sea populations, IBD did not explain the genetic structure. Results instead pointed to effects of ecotype (IBE) and spatial sorting (Berggren, Tinnert, & Forsman, 2012; Shine, Brown, & Phillips, 2011), and indicated that IBT might also influence the neutral genetic structure in pike. For adaptive differentiation, temperature and salinity appear to be two key environmental factors driving local adaptations (IBA). Below, we discuss our findings in relation to previous studies and predictions from theory, and their implications for management.

Within population diversity

Observed within population heterozygosity (HO ) for the full dataset (0.17 - 0.28, Table 1 ) were considerably lower than reported in previous microsatellite studies of Baltic Sea pike (0.22 - 0.66; Bekkevold et al., 2015; Nordahl et al., 2019; Sunde et al., 2020a; Wennerström et al., 2016). This may in part reflect the use of different markers. This interpretation is supported by the results from the study by Sunde et al. (2020a), in which both microsatellites and RADseq were utilized to study genetic diversity and structure of three pike populations (that are also included in the present study). The results revealed similar estimates of SNP heterozygosity as in the present study, and microsatellite heterozygosity estimates that are comparable to those obtained in the other studies utilizing microsatellites (0.40 - 0.57). Regardless of marker type, both the present (Table 1 ) and previous studies show that within-study estimates of HO andHE are similar, which indicates that the populations do not interbreed to a large extent (indicated byHO > HE ), and do not show any clear signs of inbreeding (indicated byHO < HE ).
Diversity is imperative for eco-evolutionary success. Theory and empirical evidence concur that larger genetic and phenotypic variation among individuals within populations may promote establishment success, stabilize population dynamics, allow for faster range expansions, and reduce extinction risk (Des Roches et al., 2018; Forsman, 2014; Forsman & Wennersten, 2016; Hughes, Inouye, Johnson, Underwood, & Vellend, 2008; Waldman, Wilson, Mather, & Snyder, 2016). Our present estimates do not indicate any clear differences in intrapopulation diversity among the ecotypes, neither for estimates based on the full nor on the adaptive dataset (Table 1 ). However, within ecotypes, some populations seem to harbor larger genetic diversity than others, which might reflect more heterogeneous environments or result from a higher degree of gene flow. Yet, these conclusions must be tentative because of the relatively low number of adaptive loci and the unbalanced sampling design, with varying numbers of populations representing the three ecotypes.

Allopatric freshwater populations were strongly neutrally differentiated

We found that all populations were significantly differentiated (based on FST ; Figure 2 ), which is consistent with previous studies (Bekkevold et al., 2015; Larsson et al., 2015; Möller et al., 2020; Nordahl et al., 2019; Sunde et al., 2020a; Wennerström et al., 2016). Levels of pairwise population differentiation ranged from low to high (FST 0.06 - 0.25), and the highest levels of differentiation were found in the comparisons including either the anadromous population Harfjärden, or at least one freshwater lake population (Figure 2 ). This, together with the finding that all three freshwater populations were strongly differentiated from each other, and distinct from the other two ecotypes point to an important effect of allopatry.
It has been suggested that pike experienced drastic population declines and/or bottlenecks following postglacial recolonization across the Northern Hemisphere (Jacobsen, Hansen, & Loeschcke, 2004), and that the succeeding isolation of local river and lake systems might have resulted in genetic drift and differentiation among freshwater populations (Bekkevold et al., 2015). Along this line of argument, the substantial differentiation among the freshwater populations likely reflects the combined effects of founder events, genetic drift, low gene flow, more distinct reproductive isolation, divergent selection, and longer time since divergence (because the freshwater ecotype is the original ecotype).

Sympatric Baltic Sea populations were also differentiated

Estimates of population differentiation in Baltic Sea pike differ markedly among studies. Notably, studies on smaller spatial scales report higher levels of differentiation (0.026 - 0.394 in Bekkevold et al., 2015; 0.013 - 0.396 in Nordahl et al., 2019) than those conducted on large spatial scales (0.005 - 0.135 in Laikre et al., 2005; and -0.003 - 0.14 in Wennerström et al., 2016). The values in the present study are relatively high compared to the other large-scale studies (0.06 - 0.22). While it is possible that the use of different markers influence the estimates ofFST , Sunde et al. (2020a) arrived at similar FST -values for microsatellites and RADseq SNPs, indicating consistency across the markers. So far, few attempts have been made to compare the resolution yielded by various genetic markers, and little is therefore known about whether the use of different marker types contributes to heterogeneity of results among studies of other organisms (Sunde et al., 2020a).
Accurate estimation of genetic differentiation requires correct population assignments of individuals. In our study system, where individuals only separate during spawning, sampling during foraging season entail the risk of erroneously grouping individuals from several populations, which results in underestimations of differentiation. In the present study, the anadromous individuals were therefore sampled in their freshwater spawning habitats to assure accurate population assignment. Previous studies differ in their sampling regimes/designs, which might have affected differentiation estimates, and complicates comparisons among studies. That we found stronger population differentiation than reported before for large-scale studies of Baltic Sea pike nevertheless indicates that the populations are more isolated than previously believed. This can have implications for management, as it indicates that each population should be considered as a separate unit. The importance of sampling scheme and associated implications for interpretation of result pertaining to population structure potentially extend also to other migrating species
(such as whitefish Coregonus maraena, Olsson, Florin, Mo, Aho, & Ryman, 2012; and perch Perca fluviatilis, Olsson, Mo, Florin, Aho, & Ryman, 2011), but to our knowledge this issue has not been systematically evaluated.
The analyses of genetic structure based on the full dataset showed signs of genetic clustering for the anadromous ecotype (Figure 1 and S4-S5 ), which is indicative of gene flow and/or recent divergence. All populations that were assigned to the shared ’anadromous genetic cluster’ spawn in localities on the Swedish mainland, and even the population from Ängerån (which resides in the north of the Baltic Sea) was assigned to this cluster. The two anadromous populations that were not assigned to the shared cluster (Askeby from Denmark and Harfjärden from the island of Öland) were the only two anadromous populations not spawning in localities on the Swedish mainland. Genetic clustering of populations along the coast have also been reported in previous large-scale studies of pike (Laikre et al., 2005; Wennerström et al., 2016). It has also been shown for pikeperch (Sander lucioperca ) in the northern part of the Baltic Sea in a study by Säisä, Salminen, Koljonen, and Ruuhijärvi (2010), who showed that coastal populations formed one genetic cluster, whilst freshwater lake populations showed strong genetic differentiation and formed distinct clusters.
All the anadromous populations were differentiated from each other (Figure 2 ), and previous studies report low levels of gene flow among anadromous pike populations (Nordahl et al., 2019; Sunde et al., 2020a; Tibblin et al., 2015). It is therefore unlikely that gene flow is sufficient to explain the clustering of anadromous populations. Coherent with previous studies based on mitochondrial DNA from pike across Northern Europe (Maes, Van Houdt, De Charleroy, & Volckaert, 2003; Skog, Vøllestad, Stenseth, Kasumyan, & Jakobsen, 2014), the results from the phylogenetic analysis revealed low levels of genetic variation and shallow branching among the populations (Figure. 5 ), which is indicative of recent divergence. However, the results did not provide any firm evidence for more recent divergence among the mainland populations. The use of more dense SNP data or longer reads might allow for the detection of clearer phylogenetic signals and higher resolution (Cariou, Duret, & Charlat, 2013).
That the anadromous populations from Öland and Denmark were distinct from the other anadromous populations may reflect different evolutionary histories, and a combination of founder events followed by divergent selection and stochastic processes. That the populations inhabiting the East Coast of Öland show strong differentiation from the mainland populations in the Kalmar sound region has been reported previously (Nordahl et al., 2019), and it has been suggested to result from the open water between the Sweish mainland and the island acting as a reproductive barrier. Similarly, the open water between Denmark and Sweden may constitute a reproductive barrier that has facilitated differentiation. In addition, previous work has shown that the population from Öland and one of the mainland anadromous populations experience different environmental conditions during spawning, and that this has resulted in the evolution of local adaptations during early fry development (temperature, Sunde et al., 2019; salinity, Sunde et al., 2018). It is therefore possible that the high level of differentiation partly reflects IBE/IBA, in addition to geographic separation. The population on Öland generally spawns earlier than the mainland populations (Sunde et al., 2019), and it is possible that IBT also has contributed to the genetic differentiation among anadromous pike populations reproducing in localities on the mainland and the island of Öland. To the extent that the timing of spawning migration is heritable (Tibblin et al. 2016), differences in timing of reproduction among populations may impair the success of individuals that attempt to spawn in a population different from where they were born, and thereby reduce gene flow. An alternative explanation is that this pattern reflects multiple evolutionary origins of anadromy. Skog et al. (2014) suggest that there are two clades of pike in the Baltic Sea, and it is possible that anadromy has evolved independently within the clades, but our present results do not provide conclusive evidence.
Taken togehter, more recent divergence of the anadromous mainland populations, different evolutionary histories, or multiple evolutionary origins all remain plausible explanations for the observed structure and differentiation, but we are unable to discriminate among them based on existing data. Future studies that include data for additional resident, anadromous and freshwater populations from other regions around the Baltic Sea are required to formally evaluate the competing hypotheses.
That the two Denmark populations (Stege Nor and Askeby) clustered together, despite belonging to different ecotypes, suggests that, in addition to ecotype, geography influence the genetic structure. Unlike previous large-scale studies of pike (Laikre et al., 2005; Wennerström et al., 2016), we found no evidence of IBD. This lack of IBD is perhaps partly explained by the long geographic distances between many of the populations. Patterns of IBD are more pronounced for more closely located populations, and gradually disappear with increasing geographic separation (Hutchison & Templeton, 1999; Meirmans, 2012; Tinnert, Hellgren, Lindberg, Koch-Schmidt, & Forsman, 2016; van Strien, Holderegger, & Van Heck, 2015). Consistent with this notion, the results from fastSTRUCTURE showed some signs of gene flow between the more closely located populations (Figure 1 ). It is therefore likely that IBD is of importance for local genetic structuring, but that other processes such as selection and drift have stronger effects on large scales.
The db-RDA for the full dataset indicated that both ecotype and latitude influence neutral genetic structuring. Results further indicate that whereas ecotype might be one of the main factors influencing neutral genetic structure, latitude appear to explain variation among populations within the anadromous ecotype.
Taken together, the findings reported in the present and previous studies suggest that the patterns of genetic structure observed in Baltic Sea pike have been shaped by an interplay between geography and divergent selection associated with the environments occupied by the different ecotypes (i.e., combined contributions of IBD, IBE, IBA, and IBT), as discussed below.

Adaptive genetic variation and structure

When the adaptive dataset (comprising outlier loci) was analysed, some contrasting patterns of structuring emerged. The clustering of the anadromous mainland populations that was evident for the neutral dataset was not present in the adaptive dataset (Figure 1 ). Instead, a main pattern of structuring associated with latitude appeared. The importance of latitude on adaptive differentiation was also evidenced by the db-RDA, which revealed a direction of separation associated with latitude that corresponded relatively well with the CAP1 axis (Figure 3b ). The effect of latitude appeared to be especially important within the anadromous ecotype, which was indicated by the direction of separation associated with latitude aligning with the separation among the anadromous individuals. This likely reflects that the anadromous populations included in this study covered a latitudinal range along the environmental clines in the Baltic Sea.
The significant interaction effect between midrange salinity and latitude indicates that the importance of these two factors differs according to the level of structuring, and further shows that salinity alone does not explain the patterns. Instead, selection associated with multiple environmental factors that co-vary with latitude (e.g.temperature and salinity) probably contributes to adaptive genetic structure. This is in agreement with previous findings that salinity and temperature regimes have resulted in locally adapted populations (Sunde et al., 2019; Sunde et al., 2018). Previous studies have also reported on evolution associated with salinity tolerance in other fish species in the Baltic Sea, including three-spined stickleback (Guo et al., 2015; Hasan et al., 2017), European flounder (Momigliano et al., 2017), and Atlantic herring (Berg, Slotte, Andersson, & Folkvord, 2019; Lamichhaney et al., 2012), emphasizing the general importance of salinity.
The role of environmental conditions in shaping adaptive structure was further supported by the finding that two of the freshwater populations (Kosta and Hamnaryd) that showed strong neutral differentiation shared an adaptive genetic cluster (Figure 1 ). This suggests that, despite being geographically and reproductively separated, similarities in environmental conditions between these two freshwater lakes have resulted in adaptive similarity via convergent evolution. The contribution of convergent evolution was further supported by the finding that the db-RDA based on the adaptive dataset showed considerably more overlap between populations (including the two freshwater populations sharing a genetic cluster in fastSTRUCTURE) and ecotypes compared to analysis of the full dataset (Figure 3 ).
The conclusion that environmental conditions, and in particular salinity, influence genetic structure was corroborated by the the outlier analyses (Table 2 and 3 ). The results specifically revealed that outliers residing in genes that have been suggested to be associated with salinity tolerance were identified in both analyses.The remainder of the SNPs might reside in genes that have not yet been identified in the pike genome, but it is also possible that these are non-coding loci linked to regions under selection.

Implications for management

Biodiversity is under threat worldwide by natural and anthropogenic environmental makeovers, climate change, and overexploitation. The level of genetic diversity within and among populations can influence the eco-evolutionary success of species, as well as the functioning of ecosystems, and this must inform management and protection of fish and the ecosystem services they provide. Genetic and phenotypic variation is required for populations to respond to selection and adapt to changing and novel environmental conditions (Charlesworth & Charlesworth, 2017; Roff, 1997; Wennersten & Forsman, 2012). There is also potential for the consequences of genetic variation to go beyond the level of the species, as it can influence community structure and ecosystem functioning (Des Roches et al., 2018; Hughes et al., 2008). Being important predators, competitors, and prey to other species, there are many ways by which pike and other species of fish can affect the functioning of lakes, rivers, coastal ecosystems and open oceans (Brodersen, Howeth, & Post, 2015; Donadi et al., 2017; Nilsson et al., 2019; Post, Palkovacs, Schielke, & Dodson, 2008; Tamario, Sunde, Petersson, Tibblin, & Forsman, 2019). There are thus several reasons as to why the genetic diversity among and within ecotypes of pike reported in this and previous studies must not be depleted.
A key challenge for conservation is to design management actions that maintain functional genetic and phenotypic diversity both within and among populations (Hutchinson, 2008; Larsson et al., 2015; Nordahl et al., 2019; Stephenson, 1999; Tamario et al., 2019; Wright, Bishop, Matthee, & von der Heyden, 2015). The rates and directions of genetic exchange between populations may be a natural outcome of dispersal or result from management actions, such as removal of migration barriers, compensatory breeding, supplementary stocking, (re-)introductions, and translocations (Gjedrem, Gjøen, & Gjerde, 1991; McClelland & Naish, 2007; McGinnity et al., 2009; Seddon, Armstrong, & Maloney, 2007). While genetic diversity is beneficial, restoration efforts may not always generate the desirable outcome (McClelland & Naish, 2007; Verhoeven, Macel, Wolfe, & Biere, 2011; Whitlock et al., 2013). Our finding in the present study that patterns of neutral and adaptive genetic diversity differed, which has also been reported in previous studies (Leinonen, O’Hara, Cano, & Merilä, 2008; Reed & Frankham, 2001), indicates that neutral variation is not necessarily reflective of adaptive variation. Given that it is adaptive, not neutral, variation that determines the evolvability of populations and influences their capability of coping with changed environmental conditions, this emphasizes the importance for management to base descisions on analyses of adaptive genetic diversity.