Population Genetic Analysis for Stock Enhancement of Silver Sea Bream ( Rhabdosargus sarba ) in Taiwan

: Stock enhancement is a method for replenishing depleted wild ﬁnﬁsh populations by supplementing them with hatchery-raised ﬁsh. In Taiwan, silver sea bream ( Rhabdosargus sarba ) is a predominant commercial species involved in stock enhancement projects. Although management agencies conduct stock enhancement projects, there are a lot of private releases without records. Stock enhancement is performed by the private aquaculture sector without accurate genetic records, potentially leading to unintended consequences for wild populations. We analyzed the genetics of 459 wild and 701 hatchery-reared specimens from nine batches produced by various hatcheries. Wild and hatchery-reared samples could be considered two separate clades by using a set of stable and informative microsatellite markers including type I (from gene introns and 3 (cid:48) UTR) and type II markers (randomly picked up from genome). Type I microsatellite markers could more sensitively reﬂect the loss of genetic diversity more than type II markers in the domestication process. All specimens were considered native by using mtDNA COI and microsatellites. The genetic composition of the wild population is relatively simple, and the estimated low contribution rate of the hatchery stock (1.3–10.9%; 6–50 / 459) indicated a weak but signiﬁcant genetic e ﬀ ect of stock enhancement. Therefore, establishing standards for the stock enhancement of silver sea bream for more e ﬀ ective supplementation of wild populations is imperative.


Introduction
Stock enhancement, the supplementation of dwindling natural populations with hatchery-reared populations, has become a frequent practice for addressing the deterioration of wild fishery resources [1][2][3][4]. In Taiwan, stock enhancement of marine species has been applied for at least 20 years through the development of hatchery-rearing technology [5,6]. However, hatchery-reared fish are usually produced from relatively few broodstock animals with limited genetic diversity. Genetic drift occurs in a hatchery population when the effective population size is small, akin to an anthropogenic founder effect [7]. Furthermore, the same broodstock is bred in the same hatcheries year after year, which makes the genetic makeup of these populations significantly differ from the natural population. The subsequent interbreeding between released and wild populations may change the genetic structure of natural Neighbor-joining (NJ) topology derived from partial COI gene of genus Rhabdosargus. The branch lengths were computed using the maximum composite likelihood method. All ambiguous positions were removed for each sequence pair (pairwise deletion option). The root was temporarily set at the deepest branch (midpoint rooting method). Numbers at nodes indicate topological support: 1000 bootstrap replicates in minimum evolution (ME), neighbor-joining (NJ) and maximum likelihood (ML) (ME/NJ/ML).  The branch lengths were computed using the maximum composite likelihood method. All ambiguous positions were removed for each sequence pair (pairwise deletion option). The root was temporarily set at the deepest branch (midpoint rooting method). Numbers at nodes indicate topological support: 1000 bootstrap replicates in minimum evolution (ME), neighbor-joining (NJ) and maximum likelihood (ML) (ME/NJ/ML).

Genetic Diversity within Populations
Across the five microsatellite markers, all the individuals were genotyped successfully. No monomorphic loci were found among the 11 populations. In total, 57 alleles were detected from the samples; the marker Pma4 exhibited the highest number of alleles per locus in all individuals (16 alleles), and the marker MSTN1 exhibited the lowest number of alleles per locus (eight alleles) in all individuals. The marker Pma4 in WP exhibited the highest number of alleles per locus (16 alleles), whereas the marker MSTN2 in HT1, HT2, HH, HR1, and HR2 showed the lowest number of alleles (two alleles) ( Table 2). Allele richness ranged from 1.928 (marker MSTN2) to 4.656 (marker Pma4) per locus in all individuals. The lowest allele richness was 1.650 in HT3 (marker MSTN2), and the highest allele richness was 6.056 in HHS (marker Pma4) ( Table 2). Mean estimates of expected heterozygosity over all loci among the 11 populations were between 0.569 and 0.647 (Table 2). WT had the highest expected heterozygosity for all the loci (H e = 0.647), whereas HT1 had the lowest expected heterozygosity (H e = 0.569) ( Table 2). Of the 55 population-locus combinations, 28 displayed deviations from the HWE significant at the p < 0.05 level ( Table 2). Strong trends of deviation were observed for specific loci (marker Cm300). No possible genotyping errors were found in all loci and only one locus (Pma4) in WP shown possible null alleles ( Table 2). All population-loci combinations displaying deviations included wild populations (WP and WT). The average F IS value among 11 populations was between −0.233 and 0.016. The highest F IS (0.016) was found in WT, and the lowest F IS (−0.233) was observed in HT1 (Table 2). MSTN1, MSTN2, and MSTN3 are the type I microsatellite markers (from introns and 3 UTR of myostatin 1 gene), and Pma4 and Cm300 are the type II microsatellite markers (randomly picked up from genome). Type II microsatellite markers are higher diversity than type I microsatellite markers, but type I microsatellite markers could reflect the loss of functional diversity sensitive more than type II markers in the domestication process (Type I: wild/hatchery = 27/14; Type II: wild/hatchery = 27/20) (Figure 2).  FIS (−0.233) was observed in HT1 (Table 2). MSTN1, MSTN2, and MSTN3 are the type I microsatellite markers (from introns and 3′UTR of myostatin 1 gene), and Pma4 and Cm300 are the type II microsatellite markers (randomly picked up from genome). Type II microsatellite markers are higher diversity than type I microsatellite markers, but type I microsatellite markers could reflect the loss of functional diversity sensitive more than type II markers in the domestication process (Type I: wild/hatchery = 27/14; Type II: wild/hatchery = 27/20) (Figure 2).

Genetic Differentiation among Populations
Pairwise comparisons between sampling localities were also performed. The pairwise F ST values ranged from 0 (HT4-HT5-HHS and WT-WP) to 0.027 (HT3-WP), suggesting a low level of genetic differentiation (Table 3). Most comparisons among hatchery populations were not significantly different, but comparisons revealed that HT2 and HT3 differed significantly (p < 0.0009) ( Table 3). WP and hatchery populations were significantly different, especially WP and HT (HT1, HT2, and HT3) ( Table 3). An initial analysis of molecular variance (AMOVA) indicated low differentiation (F ST = 0.014) with 1% genetic variation distributed among the 11 populations (Table 4). However, the overall F ST differentiation was nonetheless significant among populations, based on 999 permutations (p = 0 < 0.001; Table 4). Low genetic differentiation and high population connectivity across hatchery sampling sites were found (average F ST value = 0.008; Table 4). The F ST and gene flow (Nm) among the 11 populations were 0.014 and 17.209, respectively (Table 4).
A hierarchical AMOVA in which we defined several regions was performed, and all groupings were significantly supported by permutations (0.05 ≥ p ≥ 0.01) ( Table 5). However, when K = 2, ((HT1, Fishes 2020, 5, 19 7 of 16 HT2, HT3, HT4, HT5, HH, HHS, HR1, and HR2,) and (WP and WT)) this was determined to be the optimal grouping using the SAMOVA program and exhibited the highest intergroup variance (7.65%) ( Table 5). The best grouping result separated hatchery and wild populations. The F ST and gene flow (Nm) between hatchery and wild individuals were 0.017 and 14.829, respectively; between wild populations were 0 and infinity, respectively; and among the nine hatchery populations were 0.008 and 32.677, respectively (Table 5).  Table 5. Hierarchical analysis of molecular variance (AMOVA) of 11 populations of silver sea bream collected from hatcheries and the wild in Taiwan by using SMOVA. The best estimation of the K value (number of groups) was 2 and exhibited a stable representation (data not shown). The optimal two-cluster structure analysis was supported by the SAMOVA grouping results ((HT1, HT2, HT3, HT4, HT5, HH, HHS, HR1, and HR2) and (WP and WT)) ( Figure 3a; Table 5). All individuals from hatchery populations were grouped into the red cluster, and most individuals from wild populations were grouped into the green cluster. Only 28 individuals (6.7%) from WP were grouped into the red cluster, and these individuals were likely of hatchery origin (Figure 3a; Table 5). When clusters threshold is 70% in Structure analysis, there are six individuals grouped into the red cluster in WP (Figure 3a). Structure analysis of two groups (hatchery and wild) showed 10 hatchery individuals grouped into the green cluster (wild origin), and 41 wild individuals grouped into the red cluster (hatchery origin) (Figure 3b). We also use GenClass2 to check wild sample with Bayesian and frequency methods. There are 100 individuals grouped into the red cluster, and these individuals were likely of hatchery origin (Figure 3c). When the clusters threshold is 70% in GenClass2, there are

Discussion
Silver sea bream, Rhabdosargus sarba is a crucial aquaculture species. Aquaculture of silver sea bream began in the 1980s, and hatcheries in Taiwan are located in the Kaohsiung and Pingtung areas. Broodstocks are from the main fishery areas, which are off the Penghu Islands and west coast of Taiwan (Chiayi-Yunlin-Changhua). Stock enhancement in Taiwan began in 1973 with the building and casting of artificial reefs but no fish release. Until 1987, the stocking program was conducted with The estimated population structure based on the highest probability Structure run at K = 2. Each individual is represented by a thin vertical line, which is partitioned into K colored segments that represent individuals' estimated likelihood of membership in each of the K clusters. The red star indicates individuals with high probability of being grouped in the red cluster. When clusters threshold is 50% in Structure analysis, there are 28 individuals grouped into the red cluster and 393 individuals grouped into the green cluster in WP; when clusters threshold is 70% in Structure analysis, there are six individuals grouped into the red cluster in WP; (b) Structure analysis of hatchery and wild samples of silver sea bream in Taiwan. The estimated population structure based on the highest probability Structure run at K = 2. Each individual is represented by a thin vertical line, which is partitioned into K colored segments that represent individuals' estimated likelihood of membership in each of the K clusters. The red star indicates individuals with high probability of being grouped in the red cluster. The green star indicates individuals with high probability of being grouped in the green cluster. There are 41 individuals grouped into the red cluster in wild samples; and 10 individuals grouped into the green cluster in hatchery samples. Clusters threshold is 50% in Structure analysis. (c) There are 100 individuals from wild populations were grouped into the red cluster by using GenClass2, and these individuals were likely of hatchery origin. When the clusters threshold is 70%, there are 50 individuals grouped into the red cluster. Bayesian and frequency methods are used in GenClass2.

Discussion
Silver sea bream, Rhabdosargus sarba is a crucial aquaculture species. Aquaculture of silver sea bream began in the 1980s, and hatcheries in Taiwan are located in the Kaohsiung and Pingtung areas. Broodstocks are from the main fishery areas, which are off the Penghu Islands and west coast of Taiwan (Chiayi-Yunlin-Changhua). Stock enhancement in Taiwan began in 1973 with the building and casting of artificial reefs but no fish release. Until 1987, the stocking program was conducted with restocking broodstock (Japanese eels, Anguilla japonica) and fry (finfish) [6]. Seven finfish (5.8 million fry), four mollusk (5 million seeds), and six crustacean (30 million larvae) species were released before 1996 [6]. Records for this period are not detailed, and it was not known how many silver sea bream were released into the wild. However, we speculate that from the 1980s to 1996, there was less gene flow between the cultured and wild populations of silver sea bream, and the genetic differences between the wild and hatchery populations should have increased over six generations (2-3 years per generation). Although introducing new stock from the wild remains possible, hatcheries in Taiwan generally have maintained farmed stock since the 1980s. Many studies have demonstrated that the domestication process directly reduces genetic diversity and increases genetic differentiation among populations such as the large yellow croaker (Larimichthys crocea), barramundi (Lates calcarifer), gilthead sea bream (Sparus aurata), purple sea urchin (Paracentrotus lividus), and pearl oysters (Pinctada maxima) [25][26][27][28][29]. Porta et al. [30] noted that substantial loss of genetic variation can be found after just one generation. Therefore, although no direct evidence exists for this period, we speculate that a certain genetic difference exists between wild and hatchery populations.
Between 2002 and 2018, the official fishery organization (TFSDA) in Taiwan conducted stock enhancement of 21 species (20 finfish and one crab) with more than 133,593,000 individuals. For the official stock enhancement program for silver sea bream during 2002-2018, 5,769,700 juvenile silver sea bream were released into the waters off the Penghu Islands and west coast of Taiwan. However, because of a lack of assessment of the contribution of the release, it is unknown whether the stock enhancement has been beneficial. The TFSDA did not produce fish fry directly. Seeds for stock enhancement were provided by one or several private hatcheries, a case in which the genetic information of stock is unknown. Thus, genetic information regarding hatchery stock remains insufficient, and there is no research related to hatchery silver seabream in East Asia. We found Rhabdosargus sarba in East Asia is reciprocal monophyly and a few genetic differences between hatchery and wild populations, meaning no evidence of introduction or hybrid in Taiwan (Figures 1 and 3).
We conducted a molecular analysis of 701 silver sea bream larvae from nine batches produced by various hatcheries. We found that all hatchery-raised silver sea bream exhibited low genetic diversity (low Na) and could be considered a single cluster in structural analysis, AMOVA, and SAMOVA (Tables 1, 4 and 5; Figure 3). In the Kaohsiung hatchery populations, the genetic differences in the juveniles (HH) and broodstock (HHS), and the three batches of released juveniles (HT3, HT4, and HT5) were small (F ST = 0-0.005), indicating a high degree of homogenization. By contrast, the genetic differences between Pingtung hatchery (HT1, HT2) and Kaohsiung hatcheries (HH, HHS, HT3, HT4, and HT5) were large, and the F ST values ranged from 0.005 to 0.016 (Table 3). These indicated genetic heterogeneity of broodstock or that the released juveniles were not from the same stock. However, the average F ST among three hatchery populations was 0.008, and the gene flow (Nm) was 32.677. This result shows that hatchery species have extensive communication; thus, their genetics tend to be homogenous. The two wild populations were grouped into a single cluster in structural analysis, AMOVA, and SAMOVA (Tables 1, 4 and 5; Figure 3). The genetic differences between WP and WT were extremely low (F ST = 0), and the gene flow Nm was infinity (Figure 4), suggesting that WP and WT should be considered one population.
Silver sea bream in Taiwan has been cultured for nearly 40 years. Wild and hatchery populations can be clearly separated into two clusters with F ST values ranging from 0.010 to 0.027 (average F ST = 0.017). Although stock enhancement has been conducted for the last 15 years, the differences between farmed and wild fish remain relatively large, and the gene flow is 14.829, which is smaller than the genetic flow in cultured and wild populations (Figure 4). Among 459 wild individuals, 28 were considered to be hatchery-reared in structural analysis, whereas the remaining individuals had unequal proportions, illustrating the effects of release on the wild population. Deviations from HWE were observed mainly in the wild populations (WP and WT), indicating that one or more of the assumptions of the HWE model were not met ( Table 2). The nonrandom mating of hatchery-reared stock may result from a low-effective population size (small broodstock) or kinship mating. However, because this is usually not the case in wild populations, the HWE deviations may be explained by the Wahlund effect caused by hatchery-reared individuals that may have been used for stock enhancement or escaped. Despite at least 5,769,700 juveniles being released during 2002-2018, escapes were another potential source of continued gene flow. Cage farming aquaculture is mainly in the coastal waters of the Penghu Islands, and escapes can happen after typhoons. Farmed fish escapes may affect natural populations and interfere with the broodstock [31]. Two-to three-year-old silver sea bream can mature and undergo protandrous (male-to-female) sex changes later [18]. Escaped farmed fish are relatively big (at 1-2 years of age) and show a higher survival rate than juveniles used for stock enhancement. Moreover, escape through spawning is also possible (as with the gilthead sea bream) [32,33]. Although we cannot distinguish based on genetics fish used for stock enhancement from those that escaped hatcheries, they can all be traced to hatchery-reared stock. The only difference is that stock enhancement is active and conscious whereas escapes are unintended ( Figure 5).  Silver sea bream in Taiwan has been cultured for nearly 40 years. Wild and hatchery populations can be clearly separated into two clusters with FST values ranging from 0.010 to 0.027 (average FST = 0.017). Although stock enhancement has been conducted for the last 15 years, the differences between farmed and wild fish remain relatively large, and the gene flow is 14.829, which is smaller than the genetic flow in cultured and wild populations (Figure 4). Among 459 wild individuals, 28 were considered to be hatchery-reared in structural analysis, whereas the remaining individuals had unequal proportions, illustrating the effects of release on the wild population. Deviations from HWE were observed mainly in the wild populations (WP and WT), indicating that one or more of the assumptions of the HWE model were not met ( Table 2). The nonrandom mating of hatchery-reared stock may result from a low-effective population size (small broodstock) or kinship mating. However, because this is usually not the case in wild populations, the HWE deviations may be explained by the Wahlund effect caused by hatchery-reared individuals that may have been used for stock enhancement or escaped. Despite at least 5,769,700 juveniles being released during 2002-2018, escapes were another potential source of continued gene flow. Cage farming aquaculture is mainly Molecular markers were used to monitor and compare several genetic aspects of stock enhancement, including genetic diversity, sufficient broodstock, contribution rates, and hybridization between hatchery and wild individuals. Intraspecific genetic diversity is widely monitored using microsatellite markers derived from neutral or functional genes. In this study, we used only five microsatellite markers (two derived from neutral and three from functional genes) and obtained stable and clear results. The three microsatellite markers from functional genes did not exhibit deviations from HWE ( Table 2). Type I markers have less marker diversity, while type II markers have higher diversity. This is mainly because type II markers are neutral and non-functional, thus easily accumulating higher diversity. Although type II markers are neutral, our sample is significantly disturbed by human activities, and is reflected not in compliance with HWE. Besides, as our MSTN-derived microsatellites are found within transcribed regions of the genome these markers may be less polymorphic than those from untranscribed regions, but possible conservation of the primer sites could make them more transferable across species and reduce null alleles. As MSTN-microsatellite markers directly sample variation in transcribed regions of the genome, they may provide an estimate of functional diversity. MSTN (Myostatin) is a negative regulator of skeletal muscle growth, and is important for growth. Additionally, growth character is easy to select in breeding programs, so this functional diversity is also easily lost in the domestication process. Molecular markers were used to monitor and compare several genetic aspects of stock enhancement, including genetic diversity, sufficient broodstock, contribution rates, and hybridization between hatchery and wild individuals. Intraspecific genetic diversity is widely monitored using microsatellite markers derived from neutral or functional genes. In this study, we used only five microsatellite markers (two derived from neutral and three from functional genes) and obtained stable and clear results. The three microsatellite markers from functional genes did not exhibit deviations from HWE (Table 2). Type I markers have less marker diversity, while type II markers have higher diversity. This is mainly because type II markers are neutral and non-functional, thus easily accumulating higher diversity. Although type II markers are neutral, our sample is significantly disturbed by human activities, and is reflected not in compliance with HWE. Besides, as our MSTN-derived microsatellites are found within transcribed regions of the genome these markers may be less polymorphic than those from untranscribed regions, but possible conservation of the primer sites could make them more transferable across species and reduce null alleles. As MSTNmicrosatellite markers directly sample variation in transcribed regions of the genome, they may provide an estimate of functional diversity. MSTN (Myostatin) is a negative regulator of skeletal muscle growth, and is important for growth. Additionally, growth character is easy to select in breeding programs, so this functional diversity is also easily lost in the domestication process.
Although analysis with more markers facilitates the identification of individuals of hatchery origin, it greatly increases the cost for 1000 samples. Although a large amount of genetic data can be obtained through the use of NGS (next-generation sequencing) in recent years, this cannot avoid the overall research cost and complexity increase. Therefore, when the sample cannot be reduced, we consider that samples are more important than genetic markers. In the case of population genetic structure containing stochastic or uncertain factors (such as escape, religious release, etc.), Bayesian analysis shows the advantages. For example, when we have the sampling information of each batch Although analysis with more markers facilitates the identification of individuals of hatchery origin, it greatly increases the cost for 1000 samples. Although a large amount of genetic data can be obtained through the use of NGS (next-generation sequencing) in recent years, this cannot avoid the overall research cost and complexity increase. Therefore, when the sample cannot be reduced, we consider that samples are more important than genetic markers. In the case of population genetic structure containing stochastic or uncertain factors (such as escape, religious release, etc.), Bayesian analysis shows the advantages. For example, when we have the sampling information of each batch (hatchery) and wild samples, the accuracy of the analysis is significantly higher than only "hatchery" or "wild" information. That is, although this result may be simple now, we can improve our analysis through extra data (such as genetic markers, or hatchery information, etc.).
In stock enhancement programs, genetic analysis is possible only when hatchery stocks have genetically diverged from wild populations or when broodstock have been genotyped so that parentage-based tagging can be used to identify individuals of hatchery origin. Hatchery-reared individuals are released into the wild and are expected to be identified within the first generation, and, if possible, estimations regarding hybridization between hatchery-reared and wild individuals in the following few years should be made. This study demonstrated that hatchery-reared stocks exhibited admixture and could be grouped into one genetic clade. This clade may directly mix in hatcheries after imports of different stock or embryos from other hatcheries. In addition, before release, one batch of fry may represent a number of hatchery sources. The genetic composition of the wild population is relatively simple, and the low contribution rate of hatchery stock (6.1%; 28/459) indicated a weak but Fishes 2020, 5, 19 12 of 16 significant genetic effect after stock enhancement. The negative effects of stock enhancement are evident in many cases. The main effects include low survival rate, growth rate, reproductive fitness, and genetic diversity [34]. Few studies have provided direct evidence that wild stock has increased because of hatchery stocking [35]. Competition between wild and stocked fish may reduce wild population abundance and lead to replacement [35]. Low genetic diversity and high genetic differentiation are indicative of the risks of genetic management. These genetic issues should be carefully considered when designing a stocking plan, especially for fish from private hatcheries in Taiwan without adequate genetic information and management. In conclusion, this study proves the genetic evidence of stock enhancement in the past two decades. There are 6.1-8.9% (28-41/459) of individuals in the wild that come directly from hatchery sources or their hybrid. However, overall, we do not know whether these increase the fishery resources or just replace the composition of the wild population.

Materials and Methods
A total of 1160 specimens of silver sea bream, comprising 701 from the hatchery and 459 from the wild, were obtained from 2015 to 2016 (Table 6). Fresh specimens-at least 30 individuals for each batch, including juveniles, broodstock (S), and juveniles for releasing (R)-were sampled from three kinds of hatchery sources: (1) projected hatchery with native stock (HH and HHS); (2) private hatchery for TFSDA releasing without genetic information (HT: HT1-HT5); (3) unknown hatchery for religious releasing (HR: HR1 and HR2), and two field locations (WP: the Penghu Islands, release areas; WT: Chiayi-Yunlin-Changhua, on the west coast of Taiwan) ( Table 1; Figures 4 and 5). The geographical locations of these populations, the sampling locations with the abbreviated population names, and the sample size from each population are shown in Figures 4 and 5, and Table 1. Small pieces of muscle tissue (approximately 3-5 mm) were prepared from fresh (2% alcohol used for anesthesia) or frozen samples and transported to our laboratory for molecular studies and preserved in 95% ethanol. The standard proteinase K/phenol method was modified from an animal DNA extraction protocol. Moreover, 0.8% agarose gel electrophoresis was performed to assess DNA template quality.
The mitochondrial COI gene (Partial COI gene sequences data see the Supplementary Materials) was amplified for species identification by using primer pair 5 -CCACCGCTTAAACCTCAGC-3 and 5 -TCCGGAAGAAGCTAACAGGA-3 , yielding approximately 380 base pairs. Amplifications of the COI gene were performed using an initial denaturation at 94 • C for 3 min, followed by 35 cycles of 94 • C for 45 s, 55 • C for 30 s, and 72 • C for 45 s and a final extension at 72 • C for 3 min. DNA samples were purified using the QIAquick gel extraction kit (Qiagen, Taipei, Taiwan). Sequences were determined using the ABI 3100 (Applied Biosystems, Taipei, Taiwan) and were edited and aligned using DNA Baser (Heracle BioSoft S.R.L.). Partial COI gene sequences of genus Rhabdosargus were searched from NCBI. Minimum evolution (ME), neighbor-joining (NJ) and maximum likelihood (ML) topology were analyzed by MEGA X [36].
In this study, no useful information about microsatellite marks for silver sea bream (Rhabdosargus sarba). First, we tested 20 microsatellite marks from other sparid fish (data not shown). Then, five stable and informative microsatellite markers were used. Three were from silver sea bream MSTN1 (myostatin 1, GenBank: MT473239) gene introns, and 3 UTR: MSTN1 (intron1), MSTN2 (intron2), and MSTN3 (3 UTR). The other two were marker Pma4 from Takagi et al. [37] and marker Cm300 from Genebank (AB703235) ( Table 6). These markers can also be used in other Sparide fish e.g., Pagrus major (data not shown) [37]. Polymerase chain reaction (PCR) amplification was performed in 20-µL reaction volumes containing 5-10 ng of template DNA, 1× PCR buffer (10 mM Tris and 50 mM KCl, pH 9.0), each dNTP at 200 µM, 1.5 mM MgCl2, 0.5 U of Taq polymerase (Promega, Madison, WI, USA), and 4 pmol of each primer. Thereafter, PCR cycling was performed in an Autorisierter Thermocycler (Eppendorf, Germany) with initial denaturing at 95 • C for 2 min, followed by 30 cycles of denaturing at 95 • C for 30 s, annealing at locus-specific temperatures for 30 s, extension at 72 • C for 30 s, and a final extension at 72 • C for 10 min. Fragment analysis of PCR products was performed using  The observed genetic diversity (H o ), expected genetic diversity (H e ), and F IS were calculated by using Genalex 6.41 [38]. The chi-square Hardy-Weinberg equilibrium (HWE) test, pairwise F ST values and associated p values were performed using Genalex 6.41 [38]. To determine significance levels of F ST , multilocus genotypes were randomized between pairs of samples (9999 permutations), then significance after Bonferroni correction was calculated [39,40]. Possible genotyping errors and null alleles analyses were estimated by using Microchecker 2.2.3 [41]. To elucidate the population genetic structure from multilocus genotypes, an admixture model with correlated allele frequencies was developed using Structure v2.3.4 [42,43]. Five independent runs were performed for the entire data set for K values (numbers of groups) ranging from 1 to 9. All runs were based on 1,000,000 iterations of burn-in followed by an additional 5,000,000 iterations. The best estimation of the K value was conducted using Structure Harvester [44]. Summation and graphical representation of the Structure results were conducted using Clumpak [45]. A graphical representation of the Structure results was generated using Structure Plot v2.0 [46]. The estimated population structure was based on the highest probability Structure run at K = 2. Each individual was partitioned into two colored segments that represent individuals' estimated likelihood of membership in "wild and "hatchery" clusters. The clusters threshold was 50% and 70% in Structure analysis. Then, the wild individuals were also checked by GenClass2 with Bayesian and frequency methods. The clusters threshold was 50% and 70% in GenClass2 analysis. [47]. A hierarchical analysis of molecular variance (hierarchical AMOVA) was performed to partition the total genetic variance within and between regions as described in Excoffier et al. [48]. We used Arlequin to calculate analysis of molecular variance (AMOVA) ( Table 4). The most supported grouping can be automatically detected by using SAMOVA (based on Arlequin 3.5) (