Monitoring of genetic closed Tsaiya duck populations using novel microsatellite markers with high polymorphism.

Objective A set of microsatellite markers with high polymorphism from Tsaiya duck were used for the genetic monitoring and genetic structure analysis of Brown and White Tsaiya duck populations in Taiwan. Methods The synthetic short tandem repeated probes were used to isolate new microsatellite markers from the genomic DNA of Tsaiya ducks. Eight populations, a total of 566 samples, sourced from Ilan Branch, Livestock Research Institute were genotyped through novel and known markers. The population genetic variables were calculated using optional programs in order to describe and monitor the genetic variability and the genetic structures of these Tsaiya duck populations. Results In total 24 primer pairs, including 17 novel microsatellite loci from this study and seven previously known loci, were constructed for the detection of genetic variations in duck populations. The average values for the allele number, the effective number of alleles, the observed heterozygosity, the expected heterozygosity, and the polymorphism information content were 11.29, 5.370, 0.591, 0.746, and 0.708, respectively. The results of analysis of molecular variance and principal component analysis indicated a contracting Brown Tsaiya duck cluster and a spreading White Tsaiya duck cluster. The Brown Tsaiya ducks and the White Tsaiya ducks with Pekin ducks were just split to six clusters and three clusters when K was set equal to 6 and 3 in the Bayesian cluster analysis. The individual phylogenetic tree revealed eight taxa, and each individual was assigned to its own population. Conclusion According to our study, the 24 novel microsatellite markers exhibited a high capacity to analyze relationships of inter- and intra-population in those populations with a relatively limited degree of genetic diversity. We suggest that duck farms in Taiwan could use the new (novel) microsatellite set to monitor the genetic characteristics and structures of their Tsaiya duck populations at various intervals in order to ensure quality breeding and conservation strategies.


INTRODUCTION
The domestic duck was derived from the characteristically green-headed Mallard, Anas platyrhynchos, which is widely distributed over the northern hemisphere. Since more than 4,000 years ago, duck domestications have been occurred on many places of worldwide [1]. Asia keeps 87% of ducks of the whole world, with a tremendous variety of duck breeds, such as Shaoxing, Jingdin, Shanma, Liancheng, Bai, and Gaoyou duck in China, Tsaiya in Taiwan [2]. However, the abundant genetic diversity of ducks in Asia is under the pressure of foreign species improvement [3]. Recent, breeding trends in the production of duck in Asian countries seems to pursue higher performance in egg and meat production, and introduced such as selection, improvement and conservation of duck breeds are their major responsibilities. The relationship and characteristics of the eight LRI populations in this study are shown in Figure 1A and described as follows: The Germplasm Brown Tsaiya (GBT) line ( Figure 1E) was plotted using native Tsaiya ducks that have been collected by LRI Ilan branch since 1987. The appearances and genetic variety of this duck line have been maintained for fifteenth generations so far through the implementation of a rotational mating system [7]. LRI started a breeding project of Brown Tsaiya in 1984. The selection index for Brown Tsaiya applied at LRI aims at a maximum genetic improvement of egg number at 40 weeks without reducing egg weight or increasing body weight [6]. This high egg production line was namely Brown Tsaiya LRI 1 (BI) (Figure 1B). The second line of LRI's Brown Tsaiya, that is, the Brown Tsaiya LRI 2 (BII) ( Figure 1C), originated from the  fifth generation of BI and was selected to increase the number of fertile eggs after a single artificial insemination with pooled Muscovy semen since 1992 [8]. In 1996, LRI started to select the color of the egg shell from the eighth generation of BI and established a Brown Tsaiya line with blue egg shells which Taiwan people favor, called Brown Tsaiya LRI 3 (BIII) ( Figure  1D) [5,9]. For removed of black plumage on skin of Mule duck, the selection of White Tsaiya started in 1966. White Tsaiya were collected from duck farms on Taiwan into Ilan branch, LRI and applied in a closed herd to removed colored plumage. After 19 years, 1985, this population was named as Ilan White Tsaiya TLRI NO. 1 [6]. The White Tsaiya of Ilan branch were separated into 3 lines according to the direction selection: L101 as a germplasm preservation stock (Germplasm White Tsaiya, GWT) ( Figure 1F) was derived from original population without any selection project. L102, which was named Ilan White Tsaiya TLRI NO. 1 (WI) officially, as a dam line for crossing Muscovy drakes was derived from White Tsaiya control line and continue monitoring and selecting color of feather of mule ducks. L103 was an inbred line by full-sib mating, however this line has been discarded nowadays [6,10]. The Pekin duck (P) is the best known of all table breeds and originated in China and was introduced to the US in 1873. The Pekin duck in Taiwan were transferred from the US in 1954 as three hundred hatching eggs ( Figure  1G) [11]. The four new and two original egg-production duck breeds should be maintained and conserved for economic needs in Taiwan. Traditionally, the breeding strategies for a closed population were managed on the basis of visible phenotypic traits. However, the significance of breed is not just their specific identities. The separated genetic composition of each breed is the matter [12]. The information about genetic diversity within and between breeds and population genetic structures is very important to draw the essential outline for any appropriate conservation and sustainable management program. The increasing availability and development of molecular markers techniques is improving our ability to identify the genetic variation of breeds [13]. Six Tsaiya duck breeds, including original and selected, raised and developed in LRI, will keep pure breed raising and preserving. Among these breeds, there were just differences in color or single function. Therefore, the dissimilarity of genetic composition was supposed too limited to differentiate.
Microsatellites, also known as simple sequence repeats (SSRs), represent a class of repetitive sequences widely distributed in all genomes. They are abundant, ubiquitous, easy to automate, codominant, universal, robust, reliable and reproducible markers. Polymorphism patterns exhibited by microsatellites are greater than any other contemporary DNA molecular marker system [14]. In recent years, the polymorphism of short tandem repeats (STR) sequences is increasingly being supplemented with a different type of DNA molecular marker involving single nucleotide polymorphisms (SNPs) in DNA marker. However, SNPs panels have much lower discriminatory power than STR markers. Therefore, SNPs are a supplement for STR markers rather than a replacement [15]. It is comparative overview of SSR and SNP panels for studies of parentage and livestock animal populations of mammals, fishes and birds, the scholars pointed out that the SSR molecular marker has advantages of lower the development cost, lower running cost, lower quality of DNA required, lower quantity of DNA required, and lower technical expertise [14]. Some microsatellite markers have been developed from Tsaiya duck [16]. However, the markers with more polymorphism were needed for Tsaiya duck populations of LRI. In this study, the novelty developed microsatellite markers with high polymorphism were isolated, accompanied with previous markers [16], to analyze and monitor the genetic structure and variation of intra and inter-populations of LRI.
The fitness of Tsaiya duck to their environment is the results of natural selection for centuries in Taiwan. Before replacing them with "improved breeds", we ought to think about the irreversible loss of genetic diversity or has a certain conservation strategy to be developed in Taiwan. Before formulating these conservation strategies, it is one of the most important missions to clearly understand and define the genetic structure of the population of Tsaiya (strains) currently in Taiwan. The aim of this study was to create and collect a set of microsatellite markers with high polymorphism for the genetic monitoring and genetic structure analysis of Tsaiya (strains) duck populations.

Experiment animals and sample collection
A total of eight populations of LRI were investigated in this study. Six of these populations were Tsaiya ducks. There are two major breeds of Tsaiya duck: White Tsaiya and Brown Tsaiya duck. After a period of selection and breeding, two populations of White Tsaiya and four populations of Brown Tsaiya were raised in LRI, including the GBT ( Figure 1E Figure 1F) and Ilan White Tsaiya TLRI NO. 1, L102 (WI). Among the populations, the hybrid (H) was the progeny of GWT crossed with BI. The P ( Figure 1G) was a breed which was raised by LRI. The population was also investigated in this study as a reference. 80 GWT, 72 WI, 96 H, and 39 P individuals. For each bird, 3 mL of blood were drawn from the superficial plantar metatarsal vein or wing vein. Whole genomic DNA was then extracted with Genomic DNA Isolation Reagent (GenePure Technology CO., LTD, Taiwan) using the standard phenolchloroform method.

Polymerase chain reaction and polymorphism test
The selected highly repeated fragments were subjected to polymerase chain reaction (PCR) and polymorphism testing which would verify if the microsatellite loci could be amplified and used to show diversity in the investigated populations. Primers for loci amplification were designed using Primer-3plus [18]. CAG-tag (5'-CAGTCGGG CGTCATCA-3') or M13Reverse (5'-GGAAACAGCTATGACCAT-3') was added to the 5' end of one of each primer pair [19]. Following the protocol described by Schuelke [20], a fluorescent dyelabeled tag, as a third primer, was used with the primer pair to amplify the target fragments which were detectable on the capillary electrophoresis. Thirty unrelated Tsaiya ducks from the LRI were tested. PCR was performed on a 20 μL volume using a thermalcycler (GeneAmp PCR system 9700, Applied Biosystems, USA) containing 0.5 U Taq DNA polymerase (TAKARA, Kyoto, Japan), 1×PCR buffer (1.5 mM MgCl 2 ), 0.2 mM dNTP, 0.2 μM unlabeled primer, 0.04 μM tag-labeled primer, 0.16 μM dye-labeled tag, and 50 ng gDNA. The PCR cycling program was as follows: 95°C for 5 min, 35 cycles of 95°C for 30 s, 50°C to 65°C for 40 s, 72°C for 40 s, and a final elongation at 72°C for 7 min. The amplified microsatellite PCR products were analyzed with a DNA analyzer (ABI PRISM 3730 DNA analyzer, Applied Biosystem, USA). Allelic sizes of all loci were estimated relative to in-line GeneScan500 LIZ Size Standard marker (ABI PRISM, Applied Biosystem, USA). The fragment size was calibrated and analyzed with Peak Scanner Software version 1.0 (ABI PRISM, Applied Biosystem, USA). The loci which had an allele number (N a ) greater than two and similar annealing temperatures were selected for whole population analysis.

Statistical analysis
For each locus and population and across populations, commonly derived statistics from the microsatellite genotypic data including allele frequencies, the observed number of alleles (No), the observed heterozygosity (H O ), the expected heterozygosity (H E ), and the polymorphic information content (PIC) were calculated with the Microsatellite Toolkit. The Hardy-Weinberg equilibrium (HWE) test was performed using the GENEPOP computer program, which also was used to estimate F-statistics (F IT , F IS , and F ST ) for each locus, the pairwise F ST between populations, and the average inbreeding coefficient (F IS ). Nei's genetic distance (D A ) between populations was measured with Microsatellite Analyzer. The phylogenetic tree was generated via the PHYLIP program using the neighbor-joining (NJ) method with bootstrap test of 1,000 resamplings of loci with replacement. The genetic distances of the proportion of shared alleles (POSA) was used to estimate and draw a POSA individual phylogenetic tree [21].
A hierarchical analysis of variance was carried out to allow the partitioning of total genetic variance into components owing to region, population, and individuals. Computations were carried out using a hierarchical analysis of molecular variance (AMOVA) procedure, as implemented in the AR-LEQUIN 3.5 package. Principal component analysis (PCA) was performed with GENALEX v.6.501 software in order to spatially plot clusters and individuals based on the distance matrix with data standardization [21].
The model-based approach proposed for population structure analysis of the eight populations was carried out with the software STRUCTURE 2.3.1 [22] which assessed the genomic clustering (K) of the sample. To obtain a representative value of K for data modeling, ten independent runs were performed for each value from one to seven. The run length was set to 100,000 burn-ins following by 100,000 iterations. The ΔK estimated the most likely number of K that represented the population structure [23]. Subsequently, the CLUMPP software 1.1.2 [24] was used to demonstrate the optimal alignment of the 20 replicates for the same K value. The mean membership matrix across replicates was calculated with DISTRUCT v.1.1 [25].

Isolation of novel microsatellite Loci
Over 50 DNA fragments containing STRs were cloned. A polymorphism test was performed using the LRI populations, and 17 loci which had allele numbers more than two were picked. The 17 novel markers included three complex repeats (T1P4017, T3P1027, and T3P2090). In simple repeat markers, there was one dinucleotide unit (T1P1075), and four trinucleotide units (TYD037, TYD038, TYD042, and T1P2013), while others contained tetranucleotide units. These novel markers were then combined with the seven previously known markers [16], and together formed a set of markers for this study (Table 1).

Polymorphism, heterozygosity, and F-statistics of selected microsatellite loci
Polymorphism was clearly observed at all the microsatellite loci in all eight of the LRI duck populations. The genetic characteristics of the 24 microsatellite loci are listed in Table  2. The average number of alleles per locus (N a ) was 11.292. The actual number of alleles ranged from 2 (T3P1027) to 36 (TYD029). The average number of effective alleles per locus (N e ) ranged from 1.816 (T3P2090) to 16.588 (TYD029), with an average across loci of 5.370. The PIC value ranged from 0.375 (T3P1027) to 0.937 (TYD029), with an overall average of 0.708. All of the selected microsatellite loci in this study were sufficiently polymorphic, indicating that these loci were suitable for the genetic analysis of ducks.
The HE among the 24 microsatellite loci had a range of 0.450 (T3P2090) to 0.941 (TYD029), with the average value of H E being 0.746. The H O among the 24 microsatellite loci had a range of 0.347 (T3P2080) to 0.808 (TYD024), with the average value of H O being 0.591 (Table 2). However, there were four loci, namely, TYD021, T1P2013, T3P2090, and APT004, that significantly departed from the HWE (p<0.01).
The Wright's F-statistic values (F IS , F IT , and F ST ) for each locus are shown in Table 2. The average F IS for all the loci was 0.043, and the F IS per locus varied from -0.093 (APT020) to 0.304 (TYD025). The average F IT for all the loci was 0.233, and the F IT per locus varied from 0.095 (T1P1075) to 0.397 (TYD006). The mean F ST for all the loci was 0.197. This value implied that around 19.7% of the total genetic variation was caused by population differences and that 80.3% of total genetic variation was due to genetic differentiation among individuals within each population.

Intra-population genetic variability and Hardy-Weinberg equilibrium test
The genetic statistics relating to polymorphism, including H E , H O , PIC, the mean observed number of alleles, and the mean effective number of alleles (N e ), were calculated to estimate the allelic diversity at each locus of population. These genetic parameters across the 24 loci for the eight duck populations are listed in Table 3 (Table 2).
Among the eight populations, the P population had the highest observed mean number of alleles (MNA) (6.542), followed by the hybrid (6.500) and BIII (5.833) populations, while the WI population had the smallest observed MNA (3.708). Negative F IS values were observed only in one population (H), indicating an insufficient degree of inbreeding. The deviation from the Hardy-Weinberg proportions within populations (F IS ) varied from -0.047 to 0.103. The highest inbreeding effects were found in the BI population (0.103). The H population contained 16 loci that were significantly deviated from the HWE (p<0.01) ( Table 3), and were much higher than other populations.

Inter-population genetic variation
The values of F ST and D for each test population pair are summarized in Table 4. The F ST for each population pair was highly significant (p<0.05) The F ST values of the population pairs were varied from 0.083 (for the BI and H population pair) to 0.322 (for the WI and BII population pair). There were several pairs of population have sub-lowest F ST of which around 0.09 (BI and GBT pair, BI and BIII pair, H and BIII pair and H and GWT pair).
The genetic distances between the duck population pairs varied from 0.157 (for the BI and H population pair) to 0.582 (for the BII and WI population pair). There were relatively high genetic distance among the white duck populations (GWT, WI, and P) and between the other populations (p>0.3) except for GWT and H pair while GWT was paternal of H.

Population differentiation analysis
In this study, the results of the differentiation of the eight duck populations using the NJ tree are shown in Figure 2. In the NJ tree, the brown duck populations (GBT, BI, BII, and BIII) and the white duck populations (GWT, WI, and P) were independently grouped into two clads. The closest relationship www.ajas.info

893
Lai et al (2020) Asian-Australas J Anim Sci 33:888-901 in brown duck group was GBT and BIII. The following population was BI, but their bootstrap confidence value was only 56%. The furthermost population among brown ducks were BII. In the white duck clad, P located apart from GWT and WI, of which bootstrap value was 68%. The H populations were depicted as independent taxa.  12 185-237 57 R: AGAAAAACAGTATGAGCGAACAGGT 1) APT004, APT005, APT008, APT015, APT020, APT025, and APT031 from Hsiao et al [16].
A PCA of pair-wise genetic distances among the eight examined goose populations was used to represent the relative positions of the populations. The first (PC1), second (PC2), and third (PC3) principal components accounted for 25.95%, 23.62%, and 20.60% of the total variation, respectively ( Figure  3). The results of the PCA were similar to the phylogenetic tree drawn up via the NJ method except that GWT and WI was separated in PC3.

Population structure analysis
The results of the degree of variance in the duck populations via AMOVA are summarized in Table 5. In the analysis of all  the populations, the largest variation was found within individuals (76%), followed by variation among populations (19%). The variation among individuals within populations accounted for 5% of the total variation. In the Brown Tsaiya duck populations, the variation among individuals within the populations accounted for 13% of the total variation, while that among populations accounted for 7% of the total variation. In the White Tsaiya duck populations, the variation within individuals decreased to 70% while that among populations reached to 24%, while the variation among populations in Brown Tsaiya duck populations was only 13%.
The results of the STRUCTRE program analysis of all the populations are shown in Figure 4A. White Tsaiya and Brown Tsaiya were definitively separated into different clusters just at K = 2, while H and P became mixed groups. At K = 3, the GWT population was split into a new cluster. At K = 4, the BII population formed another new cluster. All populations were assigned into their cluster except that GBT bound to BII and H just has half BI and half GWT while K was set at six. Seven population excluding H were almost have their clusters until K = 8. No new cluster appeared, however, and the figures were not much different for values of K larger    In order to survey the structures on intra-clad basis accurately, analyses of the white ducks and brown ducks, respectively, were carried out ( Figure 4B, 4C). Among the brown ducks and P, the BII population became a subgroups first at K = 2, while the P population subsequently became a new subpopulation at K = 3. Four separately clusters were form at K = 4, that only GBT was a mixed population with the components of BI and BIII. At K = 5, the brown ducks and P were split into 5 clusters exactly. BII population tended to have two subgroups at K = 6. The white ducks were relatively simple. At K = 2, GWT formed a completed cluster and WI and P were the similar mixed clusters. The white ducks were split into 3 clusters exactly at K = 3. When K was over 4, the components of P had more complicated than GWT and WI.

DISCUSSION
In this study, we used the 24 higher polymorphism microsatellite markers (loci) to test the genetic structure of the eight duck populations. The results showed that the average of N a among known loci was 9.71±2.43, and the average of N a among novel loci was 11.94±9.20. Although the average of N a among novel loci was larger than known loci, the deviation of N a of novel loci was much large (N a was 2 to 36). Among those loci, a number of hypervariable loci with more than 20 alleles in the LRI populations, including TYD005, TYD024, and TYD029, were found. Such complex and hypervariable loci have previously been found in human CODIS loci and some canine microsatellite loci [26] as well. These complex loci might be beneficial for individual differentiation due to their excellent diversity. The P (ID) (probability of identity) values [27] of TYD005, TYD024, and TYD029 were 0.040, 0.040, and 0.048, respectively in the BI population, and their combined value was 7.68×10 -5 , which means these three markers alone could be used to differentiate more than ten thousand individuals.
To examine the diversity of the 17 novel microsatellite makers and the seven published markers in Taiwan Tsaiya ducks, several indicators were calculated from randomly selected samples from the LRI. The H O values of all the markers were higher than 0.3. A previous report suggested that microsatellite markers used in studies of genetic variation and distance should have H O values of between 0.3 and 0.8 in the population [28]. The N a , N e , and PIC values across all the loci were 11.29, 5.370, and 0.708, respectively, indicating higher genetic variability and diversity in the investigated ducks than among domestic ducks from other areas [29][30][31].   Significant deviation from the HWE was observed in only four of the 24 loci. By the STRUCTURE, the components in H population were just half of GWT and half of BI ( Figure  4A), which the parents of H, at K = 6. Also, an un-rooted individual phylogenetic tree ( Figure 5) showed that most of the individuals had been sorted to their population distinctly and that the H population was located just between the GWT and BI. These results indicated that the compilation of all the microsatellite markers used in this study were suitable for evaluating the genetic structures of Tsaiya duck populations in Taiwan.

Intra-population genetic variation and diversity
In Table 3, the results present the intra-population genetic variation in eight of the  [32]. There were all small and closed populations in LRI, however, all populations but GWT and GBT were established by selection of specific traits. GWT and GBT were derived from the original populations without any major selection since before 1980 [7]. The major causes of deficiencies of heterozygosity perhaps were selection of specific traits that accompanied specific alleles. The BII population tended to split into two Brown Tsaiya subpopulations when K = 6 ( Figure 4B). The deficiencies of heterozygosity on BII maybe combined a Wahlund effect. The first generation of BI crossed GWT have some characters predictably: the excess heterozygosity (H E = 0.676, H O = 0.707), far from inbreeding (F IS = -0.047) and departure HWE (16 of 24 markers HWE departure). In general, except for H, every Tsaiya duck population in LRI was a relatively low F IS (almost lower than 0.1), and the numbers of loci departing   from the HWE were less than 7. It can thus be concluded that the Tsaiya ducks in LRI probably have good breeding programs. The results of the intra-population diversity analyses for the Pekin populations from the LRI were similar, including relatively high variation (H O = 0.561, PIC = 0.589 and the highest observed MNA = 6.542) and F IS values of 0.08. WI was selected for crossing Muscovy drakes to improve the carcass appearances of mule ducks. WI had related low diversity (PIC = 0.461; H O = 0.482) and little more inbreeding (F IS = 0.098). This was probably caused by the multiple and longterm selection.

Inter-population genetic diversity and relationships
Brown Tsaiya and White Tsaiya ducks are the two major domestic breeds in this study, in addition of P and a hybrid population. According to the AMOVA results (Table 5), the proportion of genetic variation attributed to population differences in all populations was about 19%, while the proportion of genetic variation attributed to individuals within populations was 5% of the total genetic variation. Compared to other researches, differentiation was analyzed with 15 microsatellite loci in 26 breeds of Chinese indigenous ducks, and 12.43% of total genetic variation was attributed to differences among populations [33]; four breeds of Indian indigenous ducks were analyzed with 24 microsatellite loci, and 19.76% of the total genetic variation was attributed to differences among populations [34]. However, the Tsaiya populations of LRI were derived from unique species: Brown Tsaiya ducks. The formation of these breeds was just a series of selections and different periods of breeding. The differentiation among populations of Tsaiya ducks in LRI was larger than many different breeds in China and India. The only inference was that higher variance of markers were used in this study [35]. In addition. The variation attributed to individuals within populations was 5%. The high inter-population and low intrapopulation variations was very benefit to discriminate the different populations that even had similar appearances. The individual phylogenetic tree ( Figure 5) revealed that not only separated entirely eight populations, but the individuals were assigned exactly to their own population. Although all Brown Tsaiya populations were very close in PCA (Figure 3). Such high resolution, the 24 microsatellite loci became an excellent tool to monitor the genetic structure of each Tsaiya ducks population, and to consult a breeding program for breed conservation. When AMOVA analysis focused on Tsaiya ducks (six populations), the variance among populations (22%) were comparatively more heterogeneous than all populations (19%), and phenomenon which might have been caused by the differences between the Brown and White Tsaiya ducks. This hypothesis was confirmed with the drawing of the NJ phylogenetic tree, which showed two major clads of brown and white (including P) ducks (Figure 2), and through the STRUCTURE analysis, which revealed that the brown and white populations were divided into two parts when K = 2, with the H population split into two contributions ( Figure 4A). The PCA plot showed that the brown and white ducks were significantly separated in terms of PC1, while the hybrid ducks were located between them (Figure 3). The individual phylogenetic tree ( Figure 5) also divided the brown and white ducks into two clusters, with the hybrid ducks and P inserted between them. In our distinct analyses, the H population was always situated between two major groups or bisected individuals at K = 2 to K = 6, especially two parts just composed GWT and BI while almost populations were assigned to specific clusters. The 24 microsatellite markers used in this study have the ability to detect any crosses between two breeds.
The Brown Tsaiya populations in LRI were originated from GBT population and selected and derived to BI population, then BII and BIII populations were derived from BI. Moreover, LRI have been bred these four populations and have been kept their characters for becoming a breed. During such period, Monitoring their genetic characteristics and structures was very important. In AMOVA, the proportion of genetic variation of brown duck populations that came from population differences was only 13%, less than all Tsaiya duck populations, 22%. The four populations were close and compact as a group that far from other populations in 3D PCA. It's perhaps all BI, BII, and BIII were derived from GBT and all were selected by only one trait, moreover their separation was not for long time, the selective programs began at 1984, 1992, and 1996 respectively [6,8,9]. Such low diversity between populations and each population would breed for a breed, the monitoring and analyzing system should be high definition. The unrooted individual phylogenetic tree ( Figure  5) showed that four Brown Tsaiya duck populations were close but parted. The STRUCTURE analysis of the Brown Tsaiya ducks ( Figure 4B) also revealed that every individual that belonged to the same line was simply assigned to the same population, even when the portions in an individual were almost close to 100%. Our microsatellite markers had high definition and were suitable for analyzing the genetically close populations. Unexpectedly, the Brown Tsaiya and P populations were not isolated until K = 3. BII was isolated first at K = 2. However, P located on different clad of Brown Tsaiya ducks in NJ phylogenetic tree and separated from Brown Tsaiya duck cluster at PC2 in PCA. In our analysis, P still was a very different breed from Tsaiya duck. As the expectation, BII was split to a specific group first at K = 2, matching the results of neighbor joining phylogenetic tree and PCA that BII was a peripheral group of Brown Tsaiya ducks. It was maybe that BII raised and bred on Tainan, south of Taiwan, nevertheless, other Brown Tsaiya duck populations were on Ilan branch, north of Taiwan. The geographic sep-aration and continuous selection [9] probably were the causes of this phenomenon. As BII, BIII was selected from BI, but its variation or distance to GBT was even less than BI's. It could be caused by that genetic structure of BIII, picked with blue shell eggs, was closer to GBT. Otherwise, BIII and GBT have the highest H O (BIII = 0.631, GBT = 0.661) and PIC (BIII = 0.633, GBT = 0.623) ( Table 3). High polymorphism perhaps make both of them similar genetic content. The BII population was even split into two clusters at K = 6 ( Figure  4B). Selection in a closed population may have been the cause of this.
In contrast, belonging a clad in phylogenetic tree, the white duck populations (including P) were dispersal in PCA than brown ducks, that WI was split away from GWT at PC3 (Figure 3). The separating far away their co-original population maybe was caused by the progeny test for the white hybrid mule that needed long term selection. The diversity of WI was lost during the breeding procedure that was expressed in the lowest H O (0.482) and observed MNA (3.708) ( Table  3). Chang et al [10] reported the same trend that there were just three clearly distinguishable populations at K = 3, and that even at K = 4, two individual White Tsaiya ducks were still assigned fully to their own population ( Figure 4C). It showed that two white ducks populations were totally differentiation under analysis of our microsatellite markers. On the other hand, P had split into two sub-populations while K = 4.

CONCLUSION
The present study verified that the 24 high polymorphic microsatellite loci, including the 17 novel loci isolated in this study, were useful in studying the relationships and genetic diversities among the Tsaiya duck populations in LRI even there were almost similar genetically content. In addition, these markers could be applied continually in genetic monitoring for quality analysis and breeding of ducks in Taiwan. Analysis through 24 microsatellite markers, the intra-population with high diversity and low F IS , inter-population with low genetic distance and POSA or model-based approach with proper assigning and almost single portion individuals indicated the breeding conditions of Tsaiya ducks in LRI is generally fit demand. At last, the suggestion was that LRI and other duck farms in Taiwan could use the new microsatellite set to genetic monitor their Tsaiya populations at intervals for quality breeding and conservation strategies.

CONFLICT OF INTEREST
We certify that there is no conflict of interest with any financial organization regarding the material discussed in the manuscript.