Comparing RADseq and microsatellites for estimating genetic diversity and relatedness — Implications for brown trout conservation

Abstract The conservation and management of endangered species requires information on their genetic diversity, relatedness and population structure. The main genetic markers applied for these questions are microsatellites and single nucleotide polymorphisms (SNPs), the latter of which remain the more resource demanding approach in most cases. Here, we compare the performance of two approaches, SNPs obtained by restriction‐site‐associated DNA sequencing (RADseq) and 16 DNA microsatellite loci, for estimating genetic diversity, relatedness and genetic differentiation of three, small, geographically close wild brown trout (Salmo trutta) populations and a regionally used hatchery strain. The genetic differentiation, quantified as F ST, was similar when measured using 16 microsatellites and 4,876 SNPs. Based on both marker types, each brown trout population represented a distinct gene pool with a low level of interbreeding. Analysis of SNPs identified half‐ and full‐siblings with a higher probability than the analysis based on microsatellites, and SNPs outperformed microsatellites in estimating individual‐level multilocus heterozygosity. Overall, the results indicated that moderately polymorphic microsatellites and SNPs from RADseq agreed on estimates of population genetic structure in moderately diverged, small populations, but RADseq outperformed microsatellites for applications that required individual‐level genotype information, such as quantifying relatedness and individual‐level heterozygosity. The results can be applied to other small populations with low or moderate levels of genetic diversity.

Highly polymorphic microsatellite markers, given that they are available for the target species, can typically resolve the genetic structure of populations reliably even among closely related populations, and provide information on population genetic diversity, average kinship and effective population size (N e ). However, microsatellite markers have limitations such as risk of homoplasy for allele size, the presence of null alleles (Putman & Carbone, 2014;Zhang & Hewitt, 2003) or a potentially insufficient number of polymorphic loci in in the study species, which may limit their resolution power.
Single nucleotide polymorphisms (SNPs) obtained by restrictionsite-associated DNA sequencing (RADseq) producing thousands of loci thus provide an appealing alternative, particularly for species without prior genetic information, or for species and populations known to have limited amount of microsatellite variation because of prior population bottlenecks. Likewise, large panmictic populations, for which the differentiation levels are low, are especially challenging for genetic analysis. Consequently, the RADseq approach has rapidly gained popularity (Andrews, Good, Miller, Luikart, & Hohenlohe, 2016;Davey & Blaxter, 2010).
The advantages of SNP markers over microsatellites include their suitability for comparisons of both strongly and weakly diverged populations, and even species, and in revealing ancestral patterns of genetic structuring compared to microsatellites due to the slower mutation rate of SNPs compared to microsatellite regions (Andrews et al., 2016;Zhang & Hewitt, 2003). In addition, the RADseq approach can provide more reliable inferences on population structure (Bruneaux et al., 2013) and improved resolution for data sets with fewer individuals compared to the microsatellite approach (Jeffries et al., 2016). Likewise, Bradbury et al. (2015) demonstrated that SNPs obtained by RADseq were more accurate than microsatellites for characterizing introgression between Atlantic salmon (Salmo salar) from the East and West coasts of the Atlantic Ocean. Despite these advances, more information on wider range of species is still needed to compare the performance and cost-efficiency of these two marker types in determining population structure especially in small populations in need of conservation actions.
In addition to conservation applications focusing on populationlevel metrics, individual-based metrics, including relatedness, genetic diversity and family structure (full-sib and half-sib information) are valuable for managing hatchery breeding strategies, and understanding demographics (Hauser, Baird, Hilborn, Seeb, & Seeb, 2011;Stadele & Vigilant, 2016) and diversity-fitness-correlations (Hedrick & Kalinowski, 2000) in wild populations. While the fast mutation rate and high polymorphism of microsatellites allow for resolving fine-scale population structuring (Putman & Carbone, 2014), they may be less suitable for inferring genome-wide or individual-level patterns in genetic diversity (Väli, Einarsson, Waits, & Ellegren, 2008), as the number of sampled loci may not be sufficient to represent the total genome of an individual. It has also been demonstrated that SNPs obtained by RADseq produce more precise estimates of relatedness than microsatellites in a range of bird species (Thrasher, Butcher, Campagna, Webster, & Lovette, 2018). Microsatellites may be less efficient for identifying relatives particularly in populations with prior population bottlenecks and lack of gene flow, which limit allelic diversity. Further, the lack of commonly shared loci across species may be a more serious limitation for the use of microsatellite approach compared to RADseq analysis for, e.g., phylogenetic studies (Eaton & Ree, 2013;Near et al., 2018).
Resolving the relationship between individual-level genetic diversity (i.e., heterozygosity) and fitness is a long-standing question in conservation and evolutionary biology: lower diversity is expected to contribute to lower fitness (Hedrick & Kalinowski, 2000).
However, the genetic background of populations can influence the observed correlations (e.g., Tiira et al., 2006;Velando, Barros, & Moran, 2015), as can the type of genetic marker used (Miller et al., 2014;Väli et al. 2008). Moreover, published estimates of heterozygosity-fitness correlation (HFC) often have low correlation coefficients (Chapman, Nakagawa, Coltman, Slate, & Sheldon, 2009). In order to reveal a relevant HFC, the diversity of the applied genetic markers needs to represent individual genetic diversity across large areas of the genome, which is often not true for microsatellite panels (Fischer et al., 2017;Väli et al., 2008) as they usually only represent the most variable loci. Consequently, generally lower HFC has been found using microsatellite loci compared to the more numerous SNP loci; there was an almost five-fold increase in HFC when measured by RADseq approach in comparison to 10 microsatellite loci in an endangered species with low genetic diversity, the harbour seal (Phoca vituline; Hoffman et al., 2014). Furthermore, the minimum number of SNPs required for a reliable estimate of individual heterozygosity can vary between populations (Miller et al., 2014), but there is a lack of studies that have included both several populations and a large numerical range of loci in this evaluation (but see Fischer et al. (2017) for a pool-Seq approach).
Many salmonids provide excellent examples of systems where geographically connected (even sympatric) populations can be genetically isolated (e.g., Castric, Bonney, & Bernatchez, 2001;Estoup et al., 1998;Vähä, Erkinaro, Niemelä, & Primmer, 2007). Due to tremendous changes in their native breeding habitats, including the construction of dams, and overfishing particularly in the feeding areas, a large number of salmonid populations have become extirpated or declined dramatically (e.g., Bradford & Irvine, 2000;Morita & Yamamoto, 2002). As the stocking strategies (restoration or enhancement releases or both) of hatchery fish to maintain some of the impacted populations continue to be optimized, it remains necessary to characterize the best approaches for evaluating differences between hatchery brood stocks and native stocks. Further, the relatively high costs for DNA sequencing and library preparation in RADseq, as well as the potential challenges of obtaining numerous individuals for genotyping from small populations call for assessments whether there is cost-efficient increase in resolution to be gained by using RADseq analysis over microsatellites when using a low number of individuals.
In this study, we analyzed the genetic structure and diversity of three wild and one captive-bred population of brown trout (Salmo trutta L.) with both a RADseq approach and a DNA-microsatellite panel commonly used in brown trout population genetic research (e.g., Debes, Gross, & Vasemägi, 2017;Koljonen, Janatuinen, Saura, & Koskiniemi, 2013, Koljonen, Gross, & Koskiniemi, 2014Swatdipong, Vasemägi, Niva, Koljonen, & Primmer, 2010). Our applied goal was to support making informed management decisions locally. The main methodological goal was to compare the performance of microsatellite and SNP markers especially in the estimation of individual-level genetic diversity and relatedness to provide future reference for studies seeking to tailor the methodology to suit a given purpose. The latter was founded in the 1970-1980s.

| Sample collection
Total DNA was extracted from fin clips preserved in pure ethanol or dried scales using the Omega bio-tek E.Z.N:A Tissue DNA kit or Macherey-Nagel NucleoSpin Tissue kit. The quality of total DNA was controlled with electrophoresis on a 1% agarose gel and with fluorometric measurements using Qubit 2.0 with Qubit ® dsDNA HS Assay Kit (ThermoFisher Scientific).

| Microsatellite analysis
Allelic variation was determined at 16 microsatellite loci (Supporting information Table S1) for all 120 individuals (Table 1). For each sample, two multiplex PCR reactions were performed using the Qiagen Type-it Microsatellite kit in a 10-μl reaction volume with 3 μl of extracted DNA, 5 μl of kit master mix and primers with concentrations and dyes as presented in Supporting information. PCR reactions were carried out in PTC200 Thermal Cyclers (MJ Research), and the temperature profile of the PCR program was suggested in TA B L E 1 Sampling coordinates and summary of individuals used in the analysis with each marker. For RADseq data, N before filtering shown in parentheses

| Sequencing, genotyping and SNP calling
The sequencing libraries were prepared using samples outlined in Table 1 and Supporting information Table S2. From each sample, 100 ng of genomic DNA along with PstI-HF (5′CTGCAG 3′) and BamHI-HF (5′GGATCC 3′) restriction enzymes was used.
The protocol used was the same as in Lemopoulos, Uusi-Heikkilä, Process_Radtags function was used for demultiplexing, quality filtering (q) and cleaning (-c). Orthologous tags were assembled, catalogued and matched using denovo pipeline, in which the optimal parameters were obtained following Paris, Stevens, and Catchen (2017). Minimum coverage (-m), maximum mismatches between loci for a single individual (-M) and the maximum mismatches (-n) between loci for catalogue building were all set to 2. All other parameters were set to default. Population function was run for the SNP calling. On the first call, only loci that were present in at least 50% of all the individuals were kept while the rest of the parameters were set to defaults. Further filtering was done in R using the stackr (Gosselin & Bernatchez, 2016) and grur (Gosselin, 2017) packages. To exclude uninformative markers and samples with too much missing data, the number of populations and individuals where a locus had to be present were assessed based on the data (missing_visualisation function), after which new parameters were again passed onto STACKS' population function. Finally, only loci present in all four populations (-p) and in 60% of the individuals (-r) were retained. Based on this dataset, a total of five individuals with more than 20% of missing data (Table 1), and markers with more than 30% of missing data were discarded (using stackr). This was done in order to remove potential sequencing errors and uninformative missing data that could potentially bias the results (see stackr package guidelines). In addition, a filter for marker heterozygosity (maximum threshold: 0.5, as in Hohenlohe, Amish, Catchen, Allendorf, and Luikart (2011)) was applied to remove potential sequencing errors. Individual heterozygosity was between 0.14 and 0.22; thus, no individuals were excluded based on heterozygosity.
Markers were further filtered for minor allele frequency based on a local (0.02) and a global (0.005) threshold. This dataset was then assessed using missing_visualisation function and identity-bymissingness analysis with grur to confirm that populations did not cluster based on missing data. Further, only loci that were under Hardy-Weinberg equilibrium as defined by p-value threshold >0.05 in at least two populations were retained (HW tests made using pegas package; Paradis, 2010 in R).
Because the sex ratio of the samples from two of the studied populations was not known, the SNP dataset was checked for the presence of sex-linked markers to avoid introducing bias into the analysis (Benestan et al., 2017

| Genetic diversity and differentiation between populations
The following analyses were conducted on three datasets: microsatellite data from all individuals, microsatellite data from the same individuals as in the final RADseq data, and the SNPs from the final RADseq data.

| Multivariate analysis and bayesian clustering
The The likely number of distinct source populations and admixture between populations were additionally analyzed using STRUCTURE (Pritchard, Stephens, & Donnelly, 2000). STRUCTURE was repeated 20 times using a burn-in of 50,000 followed by 100,000 iterations for microsatellites and a burn-in of 200,000 and 200,000 iterations for the RADseq dataset, for each K value 2-7, after which the optimal K value was determined using CLUMPAK (Kopelman, Mayzel, Jakobsson, Rosenberg, & Mayrose, 2015), selecting the K where the mean ln likelihood converged. Once the optimal K (4) was determined, the 20 runs of STRUCTURE using K = 4 were combined using LargeKGreedy algorithm in CLUMPP v.1.1.2 (Jakobsson & Rosenberg, 2007), to which the input file was generated using STRUCTURE HARVESTER (Earl & Vonholdt, 2012). The output file of CLUMPP was visualized using DISTRUCT v.1.1 (Rosenberg, 2004).
The SNP and STRUCTURE analyses were conducted using the CSC -IT Center for Science Ltd clusters in Finland.

| Family structure and relatedness
Family structure within populations was assessed in the individuals with both SNP and microsatellite data available using COLONY v.

| Explorative multivariate analysis
Using DAPC, each population was separated from the others based on both microsatellites ( Figure 2a) and SNPs (Figure 2b). The hatchery stock (Ouv) and River Tuhkajoki population clustered together using microsatellite data from individuals included in RADseq (Figure 2c), thus displaying only three groups in total.

| Bayesian clustering using structure
The average likelihood of 20 independent STRUCTURE runs con-

| Genetic diversity and differentiation
In total, 3-18 alleles were found in the 16 microsatellite loci (the sum of all alleles was 143), with a mean of 4-6 alleles across loci within populations (Table 2) Table 2). Two individuals from the Vaarainjoki population shared an identical multilocus genotype based on 16 microsatellite markers, which is highly unusual, but may be explained by the low allelic diversity in most of the loci in the population. These individuals were full-sibs based on SNPs, but in the final set of loci their difference was <2%, which falls within the margin for genotyping error and might indicate that DNA was accidentally collected twice from the same individual, although we consider this unlikely. The average allelic richness was over three times higher for microsatellites than for SNPs (Table 2). The estimates of N e were overall low ( Table 2).
The confidence intervals for N e were overlapping and therefore did not indicate significant differences in the estimates across the datasets.
Based on the complete microsatellite data, the global estimate of F ST was 0.209 across populations, and was very similar to that based on SNPs (0.211). The pairwise F ST values between the three rivers and the hatchery stock were highly significant, and approximately at the same level in all three datasets (Table 3).

| Family structure and relatedness
According to COLONY output, the estimated number of full-sib families was similar using microsatellite and SNP markers on the TA B L E 2 Genetic diversity and effective population size (N e ) in all studied populations across the three datasets. Total number of alleles, mean locus-specific allelic richness (Ar mean) and the estimates of expected heterozygosity (H e ) and N e shown as measured from microsatellites in 120 individuals (A) and RADseq (B) and microsatellites (C) on the same 75 individuals.
LD-based estimates of N e are shown with 95% non-parametric jackknifed confidence intervals. Note the difference in scale for H e : the theoretical maximum is 1 for multiallelic microsatellites and 0.5 for bi-allelic SNPs SNPs. The average exclusion probabilities of full-sib family identification were clearly higher based on SNPs than microsatellites (Table 4), indicating SNPs contained more information to determine sib-ship. The differences between the markers were even more pronounced for half-sib assignment, where microsatellite data yielded four times more half-sib dyads compared to SNP data, with in total 76 based on SNPs and 314 based on microsatellites (Table 4). There was no correlation between the probabilities of the half-sib dyads that were identified with both markers (Figure 4a, Pearson r = −0.04, The mean difference in relatedness was low and stable between ca. 80-100 SNP loci (Figure 5b). In contrast, when looking at the microsatellite data (Figure 5a), the mean difference in relatedness compared between 1 and 16 loci was overall higher than that measured by SNPs (Figure 5b).

| Estimating individual heterozygosity and its accuracy
There was a positive, but relatively moderate correlation, between sMLH measured using 16 microsatellites and 4,876 SNPs (Figure 4b, Pearson r = 0.45, t = 4.29, df = 73, p < 0.001). The range of sMLH TA B L E 3 Pairwise F ST values for three wild brown trout river populations and one hatchery stock obtained using the full microsatellite dataset (A), and RADseq (B) and microsatellite data (C) on the same individuals TA B L E 4 The number of identified full-sib families with average exclusion probabilities and half-sib dyads with average probabilities from microsatellite and SNP data on the same individuals (COLONY software). A comparison of the probabilities of matching half-sib dyads is shown in Figure 4a Pohjajoki

(a) (b)
was much higher for microsatellites (Figure 4b) for which the average within-population standard deviation was 2.5 times higher than for the SNP markers. Further, analysis on subsets of SNP markers revealed that the individual genetic diversity was most precisely measured with 1,500 or more SNPs (mean Pearson r = 0.87 for 1,500 SNP subset, r = 0.90 for 2,000 SNP subset, Figure 6). There were no clear differences between populations in the number of SNPs required for a high correlation, but Tuhkajoki and Pohjajoki had higher variation in correlations (based on only 9 and 11 samples; Supporting information Figure S1).

| D ISCUSS I ON
We conducted a comprehensive evaluation of genetic divergence and diversity, and family structure in small, freshwater, postglacial populations of brown trout using both SNP and microsatellite markers. The results first showed that a moderately diverse microsatellite panel of 16 loci covering a total of 147 alleles can produce similar results as >4,800 SNPs obtained by RADseq for quantifying population divergence, and second, that the resolution of SNPs is higher compared to microsatellites in a multivariate analysis. Third, the results suggest that thousands of SNP markers are needed to reliably estimate individual-level heterozygosity.

| Performance of markers at population and individual levels
Moderate and high population divergence was equally well reflected by the two marker types, supporting work in other species.
For instance, in the Atlantic salmon, as few as nine SNP markers produced F ST values that were correlated to those measured using 14 microsatellites (Ryynänen, Tonteri, Vasemägi, & Primmer, 2007). Other studies have also reported very similar estimates of population divergence between less than a dozen microsatellite loci and >1,000 SNPs in European honeybees (Apis mellifera mellifera; Muñoz et al., 2017), Arabidopsis halleri (Fischer et al., 2017) and round whitefish (Prosopium cylindraceum; Morgan et al., 2017).
However, in Crucian carp (Carassius carassius), a greater isolation-by-distance was identified by RADseq than by microsatellites across Northern Europe (Jeffries et al., 2016). As in Jeffries et al.
(2016), we found stronger divergence between populations using the RADseq data in DAPC when comparing the two marker types F I G U R E 5 Pairwise differences in relatedness using subsets of loci from RADseq (a) or microsatellite data (b) from the same 75 Salmo trutta individuals. Pairwise relatedness between individuals was compared between each subset and the maximum number of loci used (100 SNP or 16 microsatellite loci) F I G U R E 6 Violin plots showing correlations between subsets of sMLH values in Salmo trutta according to different markers. For both microsatellites (msats) and SNPs, the correlation between two equal-size randomized subsets was calculated for 1,000 replicated sets of loci. Points showing means within each subset from the same individuals. A recent population genetics study on mud crab (Rhithropanopeus harrisii) also concluded that analyzing even a few individuals from each population with RADseq can be highly informative (Forsström, Ahmad, & Vasemägi, 2017).
Alternatives to RADseq-based approach include, for instance, reduced sets of SNP loci, which can perform nearly as well as thousands of loci in highly diverged populations (Henriques et al., 2018). In addition, increased resolution can be achieved by sequencing a large number of microsatellites first identified with a whole genome sequence scan (Bradbury et al., 2018).
Our results imply that one of the major advantages of RADseq over microsatellite analysis lies in the power to detect family structure within a population. This was shown by much higher probabilities of identified full-sib and half-sib families. The lower probabilities obtained with the microsatellite loci can be explained by the small N e in the studied populations. In addition, we found only ≤4 alleles in approximately half of the microsatellite loci within populations indicating relatively low diversity dominated by few alleles. Similar allelic richness as in our study have been described also in other brown trout populations from the wild and from hatcheries using partly the same loci (Aho, Rönn, Piironen, & Björklund, 2006;Koljonen et al., 2013;Swatdipong et al., 2010). Low within-population diversity could be a general phenomenon in endangered populations, suggesting that the results can be applied to other species.
Despite the higher probabilities compared to microsatellite markers, half-sib identification using SNPs could be confounded by genotyping errors, but little research has been done to investigate this thus far. However, a previous study comparing a custom-made SNP panel to 10 microsatellite loci for family identification in brown trout found relatively low overlap in full-sib identification between the markers, and more repeatable results when using ca. 3,800 SNPs than when using microsatellites (Linlokken, Haugen, Mathew, Johansen, & Lien, 2016). In contrast, 14 microsatellites and 1,728 SNPs agreed on 98% of full-sibs identified from 255 individuals in brown trout from River Altja (Ahmad, Debes, Palomar, & Vasemägi, 2018;Debes et al., 2017). The different results on family identification with the two marker types can therefore be partly explained by sample size; in this study, <30 individuals from each population were genotyped using both markers, and Linlokken et al. (2016) used 47-48 individuals per population from three populations. In addition, the frequency and diversity of alleles in each population finally determines the ability of microsatellites to correctly identify family structure. Notably, the use of microhaplotype markers, i.e., markers consisting of regions carrying several SNPs, has emerged as an advantageous approach for population genetic studies and particularly relatedness inference (Baetscher, Clemento, Ng, Anderson, & Garza, 2018). Such approach can be performed with RADSeq data and it has already been used in salmonids, e.g., for Chinook salmon stock identification (McKinney, Seeb, & Seeb, 2017).
Accurate quantification of genetic diversity at the individual level is crucial for HFC studies, which are comparing genetic diversity at specific loci or at the genome-wide level to variation in fitness-related traits (Chapman et al., 2009). HFC studies evaluate if the decreased genetic diversity (i.e., inbreeding depression) can reduce fitness in wild populations or make them more vulnerable to environmental disturbances due to lower adaptability to novel challenges (Hedrick & Kalinowski, 2000;Willi, Buskirk, & Hoffmann, 2006). Before SNP markers became widely available, HFC studies typically evaluated links between fitness and at heterozygosity at single or few highly variable loci (Balloux et al., 2004).  (Fischer et al., 2017).
While the most precise SNP subsets (>1,500 loci) in our study reached ~0.9 heterozygosity-heterozygosity correlations, perfectly evaluating the accuracy (i.e., a correlation of 1.0 across all loci) of the sMLH estimates would, however, require whole genome data. Overall, both our results and published work indicates that thousands of SNPs are needed to provide accurate estimates of genome-wide diversity, while estimates based on a few hundred or less SNP markers should be used cautiously.

| Implications for conservation and management of brown trout
All the rivers from which the study populations originated have undergone some alterations in the environment due to logging, dams, mining or forestry-induced decrease in water quality. These effects combined with historically high fishing pressure have further depressed the number of spawning individuals in the study populations. This explains the low genetic diversity and the N e estimates that were below 25 according to SNPs. In previous studies on other focal systems, estimates of N e below 50 have often been observed, which are in line with our observations (Linlokken et al., 2016;Linlokken, Johansen, & Wilson, 2014;Sonstebo, Borgstrom, & Heun, 2007;Vøllestad, 2017) suggesting that many brown trout populations in brooks are naturally small and strongly genetically differentiated. Each of the studied wild populations appeared to be an isolated unit with very limited gene flow from the other rivers or from the hatchery-reared fish, which is in line with their historical and present-day connectivity. Although Pohjajoki and Vaarainjoki are both connected to Lake Oulujärvi,  (Hansen, 2002;Hindar, Ryman, & Utter, 1991;Ozerov et al., 2016;Salminen, Koljonen, Säisä, & Ruuhijarvi, 2012), as is the case in Lake Oulujärvi area. Further, introgression from hatchery-reared fish can occur at sites distantly located to the stocking location (Vasemägi, Gross, Paaver, Koljonen, & Nilsson, 2005;Finnegan & Stevens 2008).
Wild brown trout populations may be locally adapted to their native environments (Jensen et al., 2008;Westley, Ward, & Fleming, 2013) and introductions of genetically differentiated hatchery-reared fish could negatively impact such adaptations (Reed et al., 2015) in case they lack the same adaptative characteristics. Additional genetic studies are required to understand the risks and potential benefits of interbreeding between the resident wild stocks and the migratory stock maintained in captivity and used for stockings and enhancement projects in the region.
Our results suggest that the existing wild populations in the region cannot be supported with the hatchery brood stock without risking the unique genetic composition in the wild populations. Thus, the management decisions on the wild brown trout populations need to balance the low diversity and N e indicating a high extinction risk with the conservation of seemingly unique genetic composition.
Particularly, large enough unique populations should not be mixed with hatchery stocking in line with the concept of Evolutionarily Significant Units proposed by Fraser and Bernatchez (2001), while the extremely small populations might benefit from the increase of genetic variation from controlled stocking with regional strains.

| Recommendations for conservation genetics studies
The resolution power of the marker depends on the number and frequency of alleles and loci available for each population.
Microsatellites can usually carry a several-fold higher number of alleles than SNPs, but the number of loci can be increased enormously in RADseq analysis compared to microsatellites. Thus, the abundance of SNP loci creates an inevitable advantage in resolution power, especially in cases where microsatellites may not function optimally due to high relatedness within populations or population bottlenecks. Consequently, although microsatellites have been central for characterizing post-glacial phylogenetic relationships in salmonids (e.g., Hansen, Mensberg, & Berg, 1999;Koskinen, Knizhin, Primmer, Schlotterer, & Weiss, 2002;Säisä et al., 2005), and many studies report comparable estimates of genetic divergence between microsatellite and SNP markers, RADseq approach could be preferred over microsatellites for detailed phylogeography studies.
Our study indicates the benefits or RADseq can also be pronounced for small populations with limited genetic diversity due to prior population bottlenecks, which is often the case for endangered populations. Thus, the choice of marker to be used in future studies depends on the goals of the study as well as available resources; compared to microsatellite analysis, RADseq still carries higher per-individual costs, although the overall price is decreasing with decreasing sequencing and reagent costs. There is a relatively high labor and time cost associated with establishing a new microsatellite panel, which RADseq circumvents by having similar preparation costs regardless of prior knowledge on markers. For projects aiming solely at describing population divergence across landscapes using a multiallelic microsatellites panel previously developed for the focal species, the benefits of switching to a RADseq approach can be limited. In contrast, many detailed questions on relatedness, demography, selection, genetic diversity or fine-scale population divergence greatly benefit from an approach based on a large number of SNPs (Andrews et al., 2016). The required number of SNPs for these questions, and therefore sequencing depth and cost, varies by species and populations, but is likely lower for relatedness analysis than for, e.g., HFC studies. In conclusion, SNPs were more informative than microsatellites for describing relatedness and had much higher resolution with low sample sizes from small and isolated populations.

ACK N OWLED G EM ENTS
The study was funded by the Academy of Finland (decision #286261, #266321) and by the Natural Resources Institute Finland (Luke). E. Taskila, P. K. Korhonen, T. Laaksonen, A. Leinonen, A. Toivonen, M.