Identification of novel haplotypes and interpretation of gene flow of mitochondrial DNA control region of Eurasian otter (Lutra lutra) for the effective conservation.

The number and distribution of Eurasian otters have declined during twentieth century due to human activity and water pollution. The global conservation status of Eurasian otter is presently 'Near Threatened (NT)' and strictly protected by being listed on the international legislation and conventions. A number of studies using the mitochondrial DNA (mtDNA) control region (CR) have been conducted in order to effectively apply conservation and reintroduction programs, especially in Europe. However, aside from Europe, there have been few studies concerning genetic diversity and phylogeny of Eurasian otters. Therefore, in this study, we sequenced partial mtDNA CR sequences (232 bp) from five South Korean Eurasian otters and analyzed 27 otters originating from parts of northeast Asia (South Korea, China, Japan and Russia (Sakhalin)), and Europe. Out of 232 bp partial mtDNA CR sequences, 13 polymorphic sites (5.6%) were identified and 4 novel mtDNA CR haplotypes (Lut16-19) were discovered from 12 Eurasian otters originating from northeast Asian region. In this study, a comprehensive analysis of genetic diversity and population structure of Eurasian otter between Europe and northeast Asia continents were conducted. Of these, different past demographic histories in Pleistocene period might have largely impacted the genetic structure of each population differently. In addition, low degree of gene flow, isolation by distance (IBD) pattern from geographically wide distanced dataset and analysis of molecular variance (AMOVA) also represented distinct genetic characteristics of Eurasian otter between Europe and northeast Asia.

The Eurasian otter (Lutra lutra L., 1758) is a semi-aquatic carnivore living as a top predator in a variety of aquatic habitats, such as rivers, lakes, lagoons, coastal wet lands, and marine shores. Otter is an important indicator reflecting the health of rivers and wetlands [9, 29, 41, 59]. Previously, this species was widely distributed across the Eurasian continent and parts of North Africa [41, 49]. However, the number and distribution of Eurasian otters have declined in recent decades and its global conservation status is presently ‘Near Threatened (NT)’ [49]. A number of factors, such as reduction of fish stocks, destruction of riparian habitat, hunting, road traffic accidents, and fish traps, have been suggested that threaten otter ecology [36, 59]. However, one of the most important factors causing the decline in the number of otters is water pollution caused by organochlorines dieldrin (HEOD), DDT, DDE, and polychlorinated biphenyls (PCBs), and heavy metals [38, 42, 49, 59]. Recently, this species became strictly protected by being listed on the international legislation and conventions, such as Appendix I of CITES (Convention on International Trade in Endangered Species of Wild Fauna and Flora), Appendix II of the Bern Convention, and Annexes II and IV of the EU Habitats and Species Directives. Many Asian nations also protect this species by listing it as an endangered species [41, 49].

For the purpose of effective application of conservation and reintroduction programs for the Eurasian otter, it is necessary to determine the genetic relationship of otters across their distribution range [16]. The mitochondrial DNA (mtDNA) control region (CR) has been widely used as a phylogeographical genetic marker in conservation studies of vertebrates [4, 16, 45, 61]. In the case of Eurasian otter, 15 haplotypes across the European continent have been identified so far through genetic analyses of the mtDNA CR. Genetic diversity of this species across the European continent is quite low and two dominant haplotypes are distributed across the European continent (Lut1 and Lut3, 69 and 20%, respectively) [2, 8, 13, 16, 17, 28, 40, 41, 45, 59]. The phylogeographic network from those studies is represented as a ‘star-like structure’, implying an extensive population bottleneck followed by postglacial recolonization by the surviving species from a single refugium [16, 59]. Although the genetic diversity of the European otter mtDNA CR is low, these previous studies might contribute to establishing conservation and reintroduction programs of this species.

The aim of this study is the identification of novel mtDNA CR haplotypes and their genetic diversity and the verification of obvious population genetic structure between European and northeast Asian (South Korea, China, Japan, and Russia (Sakhalin)) otters. Finally, we are willing to contribute to the otter conservation by discussing comprehensive factors that might have affected a variety of genetic aspects such as past demographic histories in Quaternary glaciations period, the role of biogeographic barriers in terms of preventing dispersal, low level of gene flow and isolation by distance (IBD) pattern of Eurasian otter between European and northeast Asian populations.

MATERIALS AND METHODS

Sampling

Five Eurasian otter carcasses originating from the Daegu metropolitan area and the adjacent Gyungsangbuk-do Province of South Korea were received from the Dong-in Animal Hospital, which is designated as a natural monument clinic by the Cultural Heritage Administration (Daejeon, South Korea). Additional mtDNA CR nucleotide sequences were obtained from GenBank (Table 1). In particular, of the northeast Asian CR sequences of Waku et al. (2016) [68], the Japanese sample LC050126, which was confirmed as sister clade of Lutra lutra clade, was not included in this study, because these sample had a long independent evolutionary history in Lutra (diverged about 1.27 million years ago), and it was genetically far more distant from the other Eurasian otter [68]. All animals in this study died of natural causes and found as a carcass in wild. Thus, we did not kill any animals for this study.

Table 1.

Details of the samples of Eurasian otters included in this study

Species	Locality	Sample	Sample type	Haplo-type	Accession number (GenBank)	Length (bp)	References
Eurasian Otter	S. Korea a)	KDO-1	Hair root	Lut16	-	1,091	This study
(Lutra Lutra)		KDO-2	Hair root	Lut17	-	1,084
		KDO-3	Hair root	Lut16	-	1,089
		KDO-4	Hair root	Lut16	-	1,089
		KDO-5	Hair root	Lut16	-	1,087
		EF672696	-	Lut16 c)	EF672696	16,505	[30] c)
		FJ236015	-	Lut16 c)	FJ236015	16,536	[27] c)
	China a)	LC049378	-	Lut18 c)	LC049378	16,537
		LC049952	-	Lut3	LC049952	16,537
		LC049377	-	Undetermined	LC049377	16,492	[68] c)
	Russia a)	LC049954	-	Lut16 c)	LC049954	16,536
	Japan a)	LC049955	-	Lut19 c)	LC049955	16,536
	Europe b)	AJ006175	-	Lut2	AJ006175	299	[8]
	Germany	AJ006177	-	Lut4	AJ006177	300
		AJ006178	-	Lut5	AJ006178	299
	U.K.	EU294257	-	Lut1	EU294257	299	[59]
		EU297256	-	Lut3	EU294256	299
		EU294255	-	Lut6	EU294255	299
	U.K./Wales	AM982528	-	Lut7	AM982528	364	[46]
	Ireland	FJ971618	-	Lut8	FJ971618	254	[17]
		FJ971619	-	Lut9	FJ971619	254
		FJ971620	-	Lut10	FJ971620	254
		FJ971621	-	Lut11	FJ971621	254
		FJ971622	-	Lut12	FJ971622	254
	Sweden	HQ113947	-	Lut13	HQ113947	300	[23]
	Northern	KC823048	-	Lut14	KC823048	300	[24]
	Fenno-Scandia	KC823049	-	Lut15	KC823049	300

a) ‘Asia’ population for the inter-population differentiation (Fst, Nm) and AMOVA analysis. b) ‘Europe’ population for the inter-population differentiation (Fst, Nm) and AMOVA analysis. c) The haplotype information of mtDNA CR sequences originating from those three studies has not been discovered and defined before we conducted haplotype analysis in this study, that is, we only borrowed mtDNA sequence information from those studies. Thus Lut16−19 is novel haplotypes for northeast Asian otters.

DNA extraction, Amplification, and Sequencing

Genomic DNA was extracted from hair roots of carcasses using a MagExtractor Genomic DNA Purification Kit (Toyobo, Osaka, Japan) according to the manufacturer’s recommendations.

MtDNA CR was amplified with primers mt-PCR-11-D-loop_F and mt-PCR-11-D-loop_R [44]. The amplification mixture contained 12.5 µl Thermo Scientific 2 × ReddyMix PCR Master Mix (1.5 mM MgCl2) (Thermo Scientific, Waltham, MA, U.S.A.), 1.25 µl forward primer (0.5 µM), 1.25 µl reverse primer (0.5 µM), 5 µl template DNA, and 5 µl nuclease free water in a 25 µl reaction. The cycling condition included an initial hot step at 95°C for 5 min followed by 35 cycles of denaturation at 94°C for 30 sec, annealing at 60°C for 1 min, extension at 72°C for 2 min and a final extension at 72°C for 10 min [44] using an Applied Biosystems 2720 Thermal Cycler (Applied Biosystems, Foster City, CA, U.S.A.).

Sequencing was outsourced to Macrogen Co. (Seoul, South Korea) with primers mt-seq-30F, mt-seq-30R, mt-seq-31F and mt-seq-31R [44].

Data analysis

A total of 27 mtDNA CR sequences were aligned using DAMBE 6 [73] and BioEdit 7.1.9 [22] with missing data and gaps ignored. Identification of polymorphic sites, analysis of DNA polymorphisms, gene flow parameter estimation, and haplotype identification were conducted with DnaSP 5.10.1 [35]. Understanding the evolutionary role of gene flow, or migration of individuals and the subsequent transfer of genes among populations, is pivotal to the management of endangered species [18, 47, 58, 60, 65]. So, gene flow was evaluated between two populations considering Wright’s F-statistics (Fst=1/ [1 + 4Nm] (if haploid, apply 2Nm instead of 4Nm)) (1931) [71] and the number of migrants per generation (Nm, where N is the local population size and m is the average rate of immigration in an “island” model of population structure) [25, 54] based on the sequence data information. By inverting Wright’s formula, the value of Nm can be estimated from Fst (Nm=[1/Fst−1]/4). The reason for estimating Nm is that this combination of parameters indicates the relative strength of gene flow and genetic drift. Genetic drift will result in substantial local differentiation if Nm <1 but not if Nm >1 [54].

Analysis of molecular variance (AMOVA) was conducted in order to evaluate the significance of haplotype distribution and to verify the population structure of this species with Arlequin v. 3.5.2.2 [15]. AMOVA was conducted between and within populations with 20,022 permutations for the significance test. The group for AMOVA consisted of two populations which are Asia, and Europe. The distance method used for AMOVA was pairwise difference with a gamma distribution (G) of 0.05.

A median-joining (MJ) network was constructed with NETWORK 5 [5]. As the median-joining (MJ) network initially constructed were quite complicated and had unnecessary median vectors (MVs) and links, we conducted the MP calculation using the NETWORK program to make the MJ network clear and concise. Best substitution patterns and rates were chosen as the Kimura 2-parameter model (+G). A discrete gamma distribution was used to model evolutionary rate differences among sites ([+G], parameter=0.0500). Molecular phylogeographic relation and evolutionary history were inferred by constructing the maximum-likelihood (ML) phylogenetic tree based on the Kimura 2-parameter model with MEGA 6 [62].

Isolation-by-distance (IBD) was evaluated in order to ascertain the existence and strength of a statistical relationship between genetic distance matrices and geographic distance matrices using Mantel test [6, 11, 20, 37, 56]. Mantel test were conducted with 10,000 permutations of randomization process for assessing significance using IBD (Isolation by Distance) v 1.52 program [6]. The geographic distance matrices were made up of pairwise distance (km) between every two location’s geographic coordinate originating from Google Maps (Google) in decimal degrees (DD) using Geographic Distance Matrix Generator (GDMG) program [14]. The genetic distance matrices were made up of Slatkin’s (1993) measure of similarity M̂ (=[1/Fst−1]/4) [55] except only for the analysis within Asian population due to the inability of producing M̂ from only one of all the pairwise Fst values, thus in case of the analysis within Asian population, we made the genetic distance matrices using pairwise Fst value. All pairwise genetic distance matrices were obtained by using Arlequin v. 3.5.2.2 [15]. Mantel test was conducted within population and among populations. The IBD analysis of within population was conducted for two groups: 1) within Asian population (n=12) which consisted of South Korea (n=7), China (n=3), Russia (n=1); Japan (n=1); 2) Within European population (n=15) which consisted of Germany (n=3), Sweden (n=3), U.K. (n=4), and Ireland (n=5), respectively. For the IBD analysis of among populations, We divided otter populations into two groups: 1) Asia (n=12), Germany (n=3)/ Sweden (n=3), U.K. (n=4), and Ireland (n=5); 2) Asia (n=12), Germany (n=3)/ Sweden (n=3), and Ireland (n=5) based on the previous studies [8, 17, 23, 24, 59]. Finally, reduced major axis (RMA) regression analysis was conducted to quantify the strength of the IBD pattern by calculating the slope and the intercept of genetic distance against geographic distance using the formulas provided by Sokal and Rohlf (1981) [57] using IBD (Isolation by Distance) v 1.52 program [6]. The scatter plot and linear regression function were visualized by Microsoft Office Excel 2007 (Microsoft Corp., Redmond, WA, U.S.A.).

Haplotype definition

A total of 15 haplotypes have been identified so far, Lut1, Lut3, and Lut6 [59] Lut2, Lut4 and Lut5 [8], Lut7 [46], Lut8, Lut9, Lut10, Lut11 and Lut12 [17], Lut13 [23] and Lut14 and Lut15 [24] based on the variable sites of partial (254–364 bp) MtDNA CR sequences of European otters (Lutra lutra). We newly discovered 4 novel haplotypes of this species originating from northeast Asian regions (South Korea, China, Japan and Russia (Sakhalin)) by analyzing partial (232 bp) mtDNA CR sequences with previous European sequences using DnaSP 5.10.1 [35]. We defined 4 novel haplotypes as Lut16 which contains KDO1, KDO3−5, EF672696, FJ236015 and LC049954, Lut17 which contains KDO2, Lut18 which contains LC049378 and Lut19 which contains LC049955.

RESULTS

Genetic variability and gene flow estimation

We obtained almost full-length mtDNA CR sequences from South Korean Lutra lutra, but the sequence length of European Lutra lutra was short (254–364 bp). Therefore, the length of the 27 aligned mtDNA CR sequences was 232 bp with missing data and gaps excluded. Total (G+C) content was 44.7%. The estimated transition/transversion bias ratio (R) was 16.17. A total of 13 variable sites (5.6%) were found, which comprised 6 singleton variable sites and 7 parsimony informative sites. Nucleotide diversity (π) was 0.01041 (standard deviation of π=0.0012) and the average number of nucleotide differences (k) was 2.416. Gene flow between northeast Asian population and European population showed relatively low degree (Fst=0.38872, Nm=0.79) due to the low Nm value (Nm<1) according to the Slatkin’s (1987) criterion [54]. Thus, genetic drift might result in local genetic differentiation between European and northeast Asian otter population.

Haplotype identification and phylogeographic relation

We identified four novel haplotypes of Eurasian otter mtDNA CR, with 12 sequences from northeast Asian regions (South Korea (7), China (3), Japan (1), and Russia (1)) in combination with 15 European otter haplotypes (Lut1−15) from the previous studies. A total of 17 haplotypes were identified based on 13 variable sites out of the 232 sites of each sequence. All substitutions were transitions except for only one nucleotide position, which showed a transversion (G↔T, at nucleotide position 109) (Table 2). Haplotype diversity (Hd) was 0.9202 for all 17 haplotypes.

Table 2.

Lutra lutra haplotypes found in Asia, and Europe. Lut16, Lut17, Lut18 and Lut19 are new haplotypes identified in this study

Median-joining network of 17 haplotypes showed moderate ‘star-like phylogeny’ (Fig. 1) and the ML phylogenetic tree with the highest log likelihood (−505.7107) is shown in Fig. 2. On the whole, the ML tree indicated relatively overall phylogeograpical pattern of Eurasian otter populations between Europe and northeast Asian region.

Fig. 1.

Median-joining (MJ) network of Eurasian otter haplotypes. Each circle represents a distinct haplotype. Note that Lut1/4/13 represents the same circle position for Lut1, Lut4, and Lut13 in this study (Table 2). MV: median vectors. The numbers on each line are polymorphic sites where substitutions occurred.

Fig. 2.

Maximum-likelihood (ML) tree of the partial mitochondrial control region sequences (232 bp) based on Kimura 2-parameter algorithms with a gamma distribution (G) of 0.05. Enhydra lutris was used as an outgroup. The number at each branch indicates bootstrap value (500 replications). Each haplotype is represented in parentheses next to each sample name and asterisk (*) indicates novel haplotypes in this study. Note that Chinese sample LC049377 represented the same haplotype information as Lut1, Lut4, and Lut13 (Table 2). But, it was hard to determine LC049377’s haplotype simply based on partial (232 bp) mitochondrial control regions sequences and their 13 polymorphic sites, so the haplotype information of LC049377 is described as ‘undetermined’ in the ML tree. Bar, nucleotide substitution.

Analysis of molecular variance (AMOVA)

Variance within populations was insignificant but variance between populations was significant, based on AMOVA and significance tests (Table 3).

Table 3.

AMOVA results and significance tests

Source of Variation a)	d.f. b)	Sum of squares	Variance components	Percentage of variation
Between Populations	1	8.424	0.56286 (Va)	37.97
Within Populations	25	22.983	0.91933 (Vb)	62.03
Total	26	31.407	1.48219
Fixation index Fst: 0.37975
Significant tests (20,022 permutations)
Va and Fst: P (random value >observed value)=0.00000
P (random value=observed value)=0.00000 P-value=0.00000 ± 0.00000

a) ‘Asia’ and ‘Europe’ populations were used for AMOVA analysis. b) degrees of freedom. Significance (P<0.0001) is indicated in bold.

Isolation by distance (IBD) and Reduced major axis (RMA) analysis

The results of IBD analysis were represented by calculating rm (Mantel correlation) value and P-value (significance level) of 10,000 permutations of randomizations. The P-value within Asia (n=12), within Europe (n=15) and Among Asia (n=12), Germany (n=3)/ Sweden (n=3), U.K. (n=4), and Ireland (n=5) were large, so IBD pattern of those groups were invalid. But, in terms of the Mantel test among Asia (n=12), Germany (n=3)/ Sweden (n=3) and Ireland (n=5), IBD pattern was valid considering that rm values were significantly different from zero (P<0.0001) (Table 4).

Table 4.

The results of Mantel test

The results of reduced major axis (RMA) analysis among Asia (n=12), Germany (n=3)/ Sweden (n=3) and Ireland (n=5), which represented valid IBD pattern, was represented based on the standard linear model formulas [57] (Table 5). The linear regression function implied strong relation (R2=0.739) between genetic distance and geographic distance (Fig. 3).

Table 5.

The results of reduced major axis (RMA) regression to calculate the intercept and the slope based on the standard linear model formulas provided by Sokal and Rohlf [57]

Fig. 3.

Reduced major axis (RMA) linear regression function and R2 value between geographic distance (km) and Slatkin’s (1993) [55] genetic distance (M̂) representing obvious isolation by distance (IBD) pattern (R2=0.739) among three pairwise genetic-geographic values consisted of otter populations originating from Asia (n=12), Germany (n=3)/Sweden (n=3), and Ireland (n=5). The slope and the intercept of RMA regression function were calculated using the standard linear model formulas provided by Sokal and Rohlf (1981) [57].

DISCUSSION

In this study, we discovered four novel mtDNA CR haplotypes from Asian regions (Lut16, Lut17, Lut18 and Lut19) in combination with the previous European otter haplotypes (Lut1−15) and verified population genetic structure and IBD pattern of Eurasian otter between Europe and northeast Asian region. In general, European otter populations went through a severe bottleneck due to Pleistocene glaciations and then recently underwent postglacial recolonization after last glacial maximum (LGM) period of the Pleistocene Era from a single refugial population [16, 17, 28, 41, 45, 59]. This severe bottleneck due to the Pleistocene glaciations and the following postglacial recolonizations from limited refugia in the Eurasian continent must have promoted genetic drifts and local genetic differentiations [26, 66, 72]. Due to these events, genetic diversity within European otter populations is substantially low [2, 8, 13, 16, 17, 28, 40, 41, 45, 59]. However, nucleotide diversity (π) and haplotype diversity (Hd) measured in this study (0.01041 and 0.9202, respectively) were relatively high in contrast to the results of previous studies within European species (π=0.0002−0.005, Hd=0.16−0.83) [8, 13, 16, 17, 59]. This phenomenon may be thought to be originating from distinct genetic aspects and different population structure between Europe and the northeast Asian otter populations. We suspected this contrasting genetic characteristic would be resulted from low degree of gene flow (Fst=0.38872, Nm=0.79) due to the limited dispersal ability of Eurasian otter and their isolation by distance (IBD) [52] between Europe and northeast Asia continent. So, we tried to evaluate IBD pattern between and within these two populations. We could detect a significant concordance between genetic distance (M̂) and geographic distance by Mantel test (rm=−0.8596, P<0.0001) and RMA regression (R2=0.739) between European and northeast Asian otter population. However, within each population, we could not detect any IBD pattern. In other words, IBD pattern was only detected from geographically wide distance dataset rather than geographically small distance dataset. Therefore, it is difficult for IBD pattern to be considered as the main factor for genetic differences between Eurasian otters in Europe and northeast Asia and in shaping the genetic structure observed in this study. Alternatively, we paid attention to past demographic histories in terms of possible glacial refugia and postglacial recolonizations of fauna and flora in Europe and northeast Asia in the Pleistocene Era. During the Pleistocene Era, the climatic oscillations in the Northern hemisphere were frequently and dramatically [1, 21, 33]. These extreme climatic changes were one of the main contributing factors for geographical distribution of vertebrate taxa, including their genetic structure and diversity [3, 51]. At that time, populations present in the refugia would be less affected by sudden climatic changes and they could survive under these harsh climatic conditions. Temperate zones, especially the peninsulas were the main refugia for many animal species [34, 64, 70]. In Europe, several peninsulas along the Mediterranean Sea, such as the Balkan, Iberian, and Italian have played a key role as refugia for many animal species [7, 32, 34, 39], and it is well known that European otter was also originating from these glacial refugia [16, 17, 28, 41, 45, 59]. In contrast, in East Asia, several past demographic and phylogeographic studies have suggested that East China, Korean peninsula, and Far-east Russia could have played as glacial refugia during the Pleistocene Era for many mammalian species [31, 33, 34, 51]. Independent demographic fluctuations and postglacial expansion must have occurred in the Eurasian otter populations in Europe and northeast Asian glacial refugia respectively, and eventually, they might have experienced different demographic histories and had dissimilar aspects of genetic diversity and population structure within each continent. Furthermore, biogeographic barriers such as mountains, rivers, seas, and deserts might have intensified these genetic differences between populations in Europe and northeast Asia concerning the low dispersal abilities of the two Eurasian otter populations. Taken together, the comprehensive factors mentioned above might result in genetic diversity and population structure presented in this study.

AMOVA and its significance test (20,022 randomizations) results showed 37.975% of genetic difference between European and northeast Asian otter populations (Fst=0.37975, P<0.0001) (Table 3), suggesting that otter populations between these two continents has substantial genetic difference.

Overall, Eurasian otter haplotypic information and 13 polymorphic sites are shown in Table 2. A total of 13 polymorphic sites showed regional differences in concordance with population genetic structure. Especially, both the 52nd (G) and the 267th (C) polymorphic sites together showed obvious differences of substitutions in two haplotypes, Lut16 (South Korea and Russia (Sakhalin)) and Lut17 (South Korea). Given the fact that South Korea and Sakhalin are quite close geographically and the phylogenetic results of the previous study [68], this pattern of distinct polymorphic sites would be able to be used in distinguishing the haplotype and the region where otters originated.

The haplotype networks in the previous European otter studies represented obvious ‘star-like phylogeny’, which indicates a severe bottleneck and recent proliferation of common haplotypes with a few mutations [16, 17, 28, 41, 45, 59]. Incidentally, our MJ network of Eurasian otter mtDNA CR haplotypes also represented moderate ‘star-like phylogeny’ on the whole. This haplotype network is in concordance with our ML tree which showed overall haplotype distribution according to the geographical origins of this species. Most of the haplotypes originating from northeast Asia (Lut16−19) have close relationships to Lut1 and Lut3, the common haplotypes in European otters (Fig. 1). Thus, despite the low level of gene flow and relatively high genetic drift, common European haplotypes (Lut1 and Lut3) may have been shared between Europe and northeast Asian continent with some sequence substitutions.

Generally, Eurasian otter has the ability to recover when habitat conditions are adequately restored and the use of corridors reconnecting scattered populations is available [28, 36]. So, where healthy populations survived, active conservation programs, which led some otter populations to expand and recover naturally, were designed aiming at improving habitat connectivity and sustaining natural dispersal through restored ecological corridors [29, 41, 48]. However, reintroduction programs of captive-reared or relocated wild otters are considered where natural colonization is not possible or natural connections are not available [41]. When considering reintroduction programs, we must avoid outbreeding depression, the reduction in fitness caused by the breakdown of coadapted gene complexes [50, 60, 63], and inbreeding depression which can threaten population viability due to the loss of genetic diversity [2, 10, 19, 41]. In this study, significant IBD pattern and distinct population genetic structure inferred from AMOVA between Europe and northeast Asian otter population were identified. Thus, when the reintroduction of captive-reared or wild otters is considered for the purpose of conservations of northeast Asian otter population as well as European otter population, the release of otters that originating from crossings between European and Asian species, known to bear mtDNA haplotypes of different origin, should be avoided to eliminate the risk of outbreeding depression and genetic introgression from captive stocks that underwent domestication [12, 41, 60, 69]. By the way, in some studies as to the IBD pattern of European otter populations for the conservative purpose, although IBD pattern existed across farther distances at the widest geographical scales, genetic sub-structure could exist at the same time because of restricted contemporary gene flow within close distances [10]. For example, although significant patterns of IBD at the widest geographical scales (e.g. ca. 600 km in Iberia, 600 km in France, 400 km in Germany and 2,000 km in Fennoscandia) were evidently verified through landscape genetic analysis, those regions were partially fragmented into the smaller geographical scales defined by the cryptic sub-group subdivisions (ranging from 140 to 1,600 km; ca. 260 and 410 km in the two populations of Iberia; ca. 120, 200 and 180 km in three populations of France; ca. 200 and 200 km in two populations of Germany; ca. 140, 1,600 and 1,600 km in the three populations of Fennoscandia). This finding means that local populations within small area mate randomly and are connected by restricted gene flow at distances not wider than a few hundred km [41]. Therefore, for the effective conservation of Eurasian otter in northeast Asian region including South Korea, we should try to find whether there will be cryptic genetic sub-structure and furthermore, genetic fragmentation which would be able to cause inbreeding or outbreeding depressions by any chance at the smaller geographical scales within the widest geographical scales where obvious IBD pattern was detected in this study.

Although this study identified 4 novel Asian otter haplotypes, and analyzed genetic diversity, population genetic structure, and IBD pattern of Eurasian otter at the widest geographical scale, there will be increasing need of studying numerous historical and current factors which will be able to influence on the genetic aspects of this species by utilizing various genetic markers in addition to the mitochondrial DNA, such as microsatellite or allozyme with extensive sampling. For the historical factors, past demographic studies such as identification of possible glacial refugia and postglacial expansions [31, 33, 34, 51], colony forming patterns [53, 66], and habitat fragmentations [10, 19, 43, 67] resulted from Quaternary glaciations across Eurasian regions, and for the current factors, studies such as different dispersal patterns during and after the range expansions [26] are additionally needed for the successful future conservation of Eurasian otter.