Population Genetic Structures of the Jellyfish Aurelia coerulea Polyps along Korean Coasts and Implications as Revealed by Mitochondrial COI

© 2021 Academia Sinica, Taiwan Open Access Population Genetic Structures of the Jellyfish Aurelia coerulea Polyps along Korean Coasts and Implications as Revealed by Mitochondrial COI Yoseph Seo1,§ , Buhari Lawan Muhammad1,§ , Jinho Chae2, and Jang-Seu Ki1,* 1Department of Biotechnology, Sangmyung University, Seoul 03016, South Korea. E-mail: [email protected] (Seo); [email protected] (Muhammad); [email protected] (Chae) 2Marine Environmental Research and Information Laboratory, Gunpo 15850, South Korea. *Correspondence: E-mail: [email protected] (Ki). Tel: +82-2-2287-5449. Fax: +82-2-2287-0070. §YS and BLM contributed equally to this work. Received 1 April 2021 / Accepted 20 July 2021 / Published 12 November 2021 Communicated by James D. Reimer The moon jellyfish, Aurelia coerulea, is is globally distributed, and its blooms have been responsible for severe environmental impacts. Benthic polyp populations are important for forming and maintaining medusa populations; however, their genetic structures are mostly unknown. Here, we analysed the genetic structure and phylogeographic pattern of A. coerulea polyps using the mitochondrial COI of 229 specimens collected from four different coastal regions of Korea. Molecular discrimination by COI assigned all polyps to A. coruelea. Population genetics revealed 53 haplotypes with high diversity and significant genetic structure, distinguishing two haplogroups (A and B) that coexist in all regions. Haplogroup A exhibited a star-like haplotype network pattern, while haplogroup B demonstrated a branched haplotype network pattern. Our results suggest that, the two haplogroups detected have existed in sympatry along Korean coasts. However, haplogroup A may have been established by a recent population expansion, while haplogroup B may have been established a long time ago. The strong genetic structure found within the polyp population of A. coerulea may have an effect on the moon jellyfish blooms on Korean coasts. Key words: A. coerulea, Benthic polyp, Cytochrome c oxidase subunit I, Genetic diversity, Jellyfish bloom. BACKGROUND Mass outbreaks of jellyfish have increased in the last several decades in near-shore regions and coastal ecosystems across the globe (Condon et al. 2012). The main causes of these increase in these jellyfish blooms are thought to be the rising of ocean temperatures (Condon et al. 2012), eutrophication (Arai 2008; Purcell et al. 2007), and the construction of marine structures (Duarte et al. 2013). The large aggregations (blooms) of jellyfish are associated with several negative socioeconomic impacts in many harbours and coastal waters around the world (Mills 2001). They have caused significant ecological and economic problems in the San Francisco Bay of Northern California (Greenberg et al. 1996), the Gulf of Mexico (Frolova and Miglietta 2020), Japan (Yasuda 2004), and South Korea (Ki et al. 2008; Yoon et al. 2018). These problems have ranged from fishery damage (Kim et al. 2012; Quiñones et al. 2013) to the clogging of seawater intakes of power plants and ships (Abdul Azis et al. 2000; Ki et al. 2008). In addition, jellyfish have also clogged the pumps used to flush the Foster City Lagoon, an artificial recreational pool, and fouled beaches (Greenberg et al. 1996). The moon jelly Aurelia is one of the predominant bloom-forming jellyfish in coastal and shelf ecosystems. Citation: Seo Y, Muhammad BL, Chae J, Ki JS. 2021. Population genetic structures of the jellyfish Aurelia coerulea polyps along Korean coasts and implications as revealed by Mitochondrial COI. Zool Stud 60:63. doi:10.6620/ZS.2021.60-63. Zoological Studies 60:63 (2021) doi:10.6620/ZS.2021.60-63 1 © 2021 Academia Sinica, Taiwan Within the genus, Aurelia coerulea von Lendenfeld 1884 is a taxon with the broadest world distribution (Dawson et al. 2005; Scorrano et al. 2017). Generally, A. coerulea is believed to have spread across oceanic regions from a single origin (Dawson et al. 2005; Ki et al. 2008). Thus, severe blooms of A. coerulea have been frequently observed worldwide (Fofonoff et al. 2018). On Korean coasts, the occurrence of A. coerulea blooms dates back to 1980. Since then, blooms have been frequently reported in harbour and port areas, especially in the western and southern coastal areas, and caused huge economic losses estimated at US$ 277.9 million per year (Kim et al. 2012). This has gained government attention, which has led to the establishment of a management policy that includes survey and removal of medusae and polyps (if necessary) of A. coerulea in Korean coastal waters (Yoon et al. 2018). However, despite enormous damage, current knowledge of the population’s genetic structures, distribution, and source of the outbreak is lacking. Without knowledge of the polyp genetic structures, geographical ranges, and their possible link to medusae blooms, it is difficult to formulate effective bloom management strategies (Dawson et al. 2015; van Walraven et al. 2016). Most studies on jellyfish blooms have focused on adult medusae. However, for many Aurelia spp., the magnitude and timing of medusae blooms are dependent on the dynamics of their benthic polyp population and the environmental factors that control it (Marques et al. 2019). In addition, the potential asexual reproductive ability of polyps is assumed to have a significant effect on the blooms of medusae (Arai 2008; Schiariti et al. 2014; Takao et al. 2014). During strobilation events, one polyp can produce ~100 ephyrae (Lucas et al. 2012); therefore, depending on the duration of the polyp stage, the asexual reproduction of the polyp can result in apparently irregular and unpredictable blooms of the medusae (Boero et al. 2008). However, although benthic populations of Aurelia spp. create large blooms of freeliving medusae, little has been investigated so far, probably because of their small body size and cryptic lifestyle, which make their detection difficult (Bayha and Graham 2009). Recent studies have focused on the benthic populations; the causes and consequences of benthic population dynamics, bloom events, and environmental tolerance have also received much attention (Lucas et al. 2012; Marques et al. 2019; Frolova and Miglietta 2020). However, studies on the population genetic structure, diversity, and phylogeography of the polyp population in Aurelia spp. are limited. van Walraven et al. (2016) were the first to use the sessile benthic stage of A. aurita to show high genetic diversity and differentiation in polyp populations similar to those in medusa populations. However, since A. coerulea is the taxon with the broadest world distribution, causing greater socio-economic damage than any other Aurelia species, it is necessary to investigate its polyp population. This will help to determine the population differentiations that may exists in the benthic stage, and their potential implications for the medusa population. Understanding the genetic structure and phylogeography of the polyp population is necessary to predict when and where jellyfish blooms will occur (van Walraven et al. 2016). The specific impacts of A. coerulea are poorly understood because it is easy to confuse this species with native cryptic species (Dawson and Jackobs 2001; Fofonoff et al. 2018). Species crypsis and morphological plasticity make the species difficult to distinguish (Dawson 2005; Scorrano et al. 2017). To resolve this, molecular tools have been applied to differentiate cryptogenic species, identify population origins, and assess the impacts of invasions. Among the molecular markers, mitochondrial cytochrome c oxidase subunit I (COI) has been widely used as a powerful molecular marker to determine the genetic variations and structure of populations in many organisms (Xu et al. 2011; Dong et al. 2015; Palraju et al. 2018; Choi et al. 2020). It is a protein-coding gene in the mitochondrial genome that serves as a useful genetic marker, providing valuable information for uncovering phylogeographic patterns and ecological adaptations (Ramšak et al. 2012). COI has been widely used for population genetic studies, especially intra-specific analysis, due to its fast evolution, high polymorphism, and easy amplification and sequencing (Hu et al. 2008; Xu et al. 2011; Palraju et al. 2018). In the present study, we used the mitochondrial COI gene to study the population genetic structure and phylogeography of A. coerulea polyp populations. The samples were collected from four different regions that cover most of the areas inhabited by A. coerulea on the Korean coastal waters. The aim of this study was to explore the genetic structure and phylogeographic patterns of the benthic populations of bloom-forming jellyfish in Korean coastal regions. This allows us to resolve the distributions and patterns of genetic differentiation among the populations of A. coerulea and identify the importance of its benthic population’s dynamics in medusae blooms. MATERIALS AND METHODS Study site and sample collection A total of 229 polyp samples of Aurelia coerulea were collected in April 2017 from four different page 2 of 13Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan geographical regions (R1, R2, R3, and R4) located on the western and southern coasts of Korea (Fig. 1). Within R1 (Jeongok, JG; Dangjin, DJ; Daesan, DS) and R2 (Baeksajang, BJ; Daecheon, DC; Galyeongdo, GL), three locations were chosen for sampling. Within R3, two locations were chosen (Usuyeong, UY; Nokdong, ND), and in R4, only one location (Jaran Bay, GS) was chosen for the sampling. Detailed information regarding Fig. 1. Sampling locations of the 229 polyp samples of Aurelia coerulea collected within four regions (R1, R2, R3, and R4) from western and southern coasts of Korea. Detailed information on the sampling locations are listed table 1. 126°18'49.0" E 128°14’05.7" E 36°54'36.8" N 50 km 0 35°17’44.9" N JG DJ DS BJ DC GL UY ND GS R1 R2 R3 R4 N Table 1. Details of the sampling locations of 229 Aurelia coerulea polyps collected from western and southern coasts of Korea Region LC N Samping date Location GPS R1 JG 15 2017.4.6 Jeongok port 37°11'20.7"N, 126°38'50.4"E DJ 19 2017.4.13 Dangjin port 36°59'12.2"N, 126°44'50.1"E DS 45 2017.5.2 Daesan port 37°00'58.7"N, 126°25'25.5"E R2 BJ 17 2017.4.13 Baeksajang port 36°35'10.9"N, 126°18'56.9"E DC 20 2017.4.14 Daecheon port 36°19'41.7"N, 126°30'32.2"E GL 23 2017.4.4 Galyeogdo port 35°43'35.9"N, 126°31'46.3"E R3 UY 3 2017.3.29 Usuyeong port 34°34'07.9"N, 126°18'32.3"E ND 33 2017.4.9~19 Nokdong port 34°31'25.8"N, 127°08'06.9"E R4 GS 54 2017.3.30 Jaran Bay 34°54'59.6"N, 128°15'11.5"E LC, location code; N, Number of sample. page 3 of 13Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan the regions and sampling locations are listed in table 1. Polyp samples were collected by scuba divers who were trained to dive for polyp sampling and distinguish polyps in situ; divers also sketched underwater polyp distribution and took photographs. Collection proceeded by scraping the polyps and their substrate using a filling knife, and then placing collected scrapings into a 1-L Kautex jar. Collected samples were immediately fixed with 100% ethanol and transported to the laboratory. Aurelia coerulea polyps were identified morphologically, and individuals were separated using tweezers under a dissecting microscope. DNA extraction, amplification, and sequencing Ethanol-preserved polyps were washed individually several times with distilled water, then genomic DNA (gDNA) was extracted using the modified cetyltrimethylammonium bromide (CTAB) method (Ausubel et al. 1989). Polymerase chain reaction (PCR) amplification of COI fragments was done by using the primer pair Fco04 5'-CAA CAA ACC ATA AAG ATA TAG GAA C-3' (Folmer et al. 1994; slightly modified) and Rco800 5'-TAT TGC TAT CAT AGC ATA AAC CAT-3' (newly designed in this study). In brief, 20 µL reaction solution was prepared for PCR. The reaction solution consisted of 2 µL template DNA, 2 µL Ex Taq buffer, 2 µL dNTP mix, 1 µL forward primer, 1 µL reverse primer, 0.2 µL Ex Taq Polymerase (TaKaRa Bio Inc., Shiga, Japan), and 11.8 µL distilled water. PCR reaction conditions were as follows: 94°C for 3 min, followed by 40 cycles of 94°C for 30 s, 55°C for 30 s, and 72°C for 1 min), and finally 72°C for 5 min. PCR results were confirmed on a 1% agarose gel using 60 µL/400 mL gel MIDORIGreen (Nippon Genetics Europe GmbH, Germany) as a fluorescent source. Confirmed samples were purified using the PCR Cleanup S & V Kit (Bionics, Daejeon, Korea) according to the manufacturer’s instructions. Purified PCR products were sequenced by Bionics Inc. (Seoul, Korea) using the same primer pair as used for PCR amplification. Sequencing was performed using the BigDye™ Terminator v3.1 Cycle Sequencing Kit (Thermo Fisher Scientific) and the synthesised sequences were analysed using an Applied Biosystems 3730xl DNA Analyser (Applied Biosystems, Foster City, CA). Sequencing data were compared with reference sequences available in the GenBank database. All sequences were identified as A. coerulea and compiled in MEGA X software (Kumar et al. 2018). COI sequences were translated into amino acid sequences using the coelenterate mitochondrial code in MEGA to confirm open reading frames. The consensus sequences were trimmed to 759-bp fragments and aligned using Clustal Omega (Sievers et al. 2011). All sequences determined in this study were deposited in the GenBank database under the accession numbers MT251296-MT251346 and MT361137-MT361314. A list of all sequences used in this study is provided in table S1. Population and phylogenetic analyses The assembled sequences were imported into Arlequin v.3.5.2.2 (Excoffier and Lischer 2010) to calculate the genetic diversity indices (number of haplotypes (Nh), number of polymorphic sites (Ps), haplotype diversity (h), and nucleotide diversity (π). In addition, Tajima’s (D) (Tajima 1989) and Fu’s (FS) scores (Fu 1997) were calculated using DnaSP v.6 (Rozas et al. 2017). To explore the phylogeographic patterns of the present A. coerulea haplotypes, we constructed a TCS haplotype network (Clement et al. 2002) for the 229 COI sequences using PopART v1.7 (Leigh and Bryant 2015). Haplotype networks were constructed using the default settings. Genetic variation among the different clusters was tested. A total of 229 COI sequences (759 bp) were imported into Arlequin v.3.5.2.2, where Analysis of Molecular Variance (AMOVA) was assessed, and a pairwise FST matrix with 10,000 permutations was constructed. The haplotype distance matrix was graphically rendered using R-3.0 software (R Core Team 2014). Principal Component Analysis (PCA) was used to further evaluate and visualise the geographic genetic structure among the populations, and was calculated in Palaeontological Statistics Software Package (PAST) v.3.25 (Hammer et al. 2001), using the haplotype frequency of the 53 unique haplotypes (Table S2). For phylogenetic analysis, we constructed a dataset of 53 COI reference sequences of the haplotypes determined in this study and 115 COI sequences of A. coerulea (Aurelia sp.1) retrieved from the GenBank database. (Table S3). Overall, 168 COI sequences (584 bp) were aligned using Sequencher v.5.1, (Gene Codes Corporation, Ann Arbor, MI), and the alignment results were trimmed to remove ambiguous regions. A phylogenetic tree was constructed using MEGA X software (Kumar et al. 2018). The tree was inferred using neighbor-joining (NJ) computed with the Kimura 2-parameter model (Kimura 1980), and taxa were clustered together with 1000 bootstrap values. The phylogenetic tree was visualised with Mega X Tree Explorer and edited with Adobe Illustrator CS6 (Adobe Systems, San Jose, CA). page 4 of 13Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan RESULTS Genetic diversity In this study, a total of 229 Aurelia coerulea specimens were analyzed. Our BLAST search comparisons for each of the 229 COI sequences revealed that our sequences shared 100% of nucleotide identity with a known COI DNA sequence from A. coerulea (GenBank accession no. LK022703, formerly Aurelia sp.1). Our genetic diversity analysis of 229 Aurelia coerulea polyps identified 51 polymorphic sites (Ps), and represents 53 haplotypes (Table 2). The overall haplotype diversity (h) and nucleotide diversity (π) across the four geographic regions were estimated to be 0.8788 ± 0.0169 and 0.005532 ± 0.003052, respectively. The highest haplotype diversity was calculated for R3 (0.9), while the lowest was for R1 (0.8079). The highest nucleotide diversity was also observed in R3 (0.0056471), and the lowest was observed in R2 (0.004577). Neutrality tests of the COI haplotypes were performed for all populations. Tajima’s D and Fu’s FS values were negative for all populations, but only Fu’s FS values at R2 and R4 were statistically significant. However, the overall values of both Tajima’s D and Fu’s FS showed significantly negative values (D = -1.484 (p < 0.1) and FS = -25.442 (p < 0.05) (Table 2). Genetic differentiation AMOVA analysis was carried out for the populations based on their geographical distance. The analysis showed that most of the genetic differentiation occurred within populations (83.87%, p < 0.001), while 18.88% of the genetic differentiation occurred among populations within regions, and -2.76% occurred among regions (Table 3). The overall FST value (FST = 0.162; p = 0) suggested moderate but significant genetic differentiation within populations with significant gene flow. Pairwise FST analysis showed the pattern of genetic differentiation for all locational populations, and analysis showed that most of the significant genetic differentiation was between location DS in R1 and the rest of the locations (Fig. 2). The largest level of genetic differentiation was found between R1 (DS) and R2 (GL), while the smallest genetic differentiation was found between location UY in R3 and other locations, except for DS in R1. However, the FST test showed no significant genetic differentiation among the rest of the locations. Phylogeography and distribution of A. coerulea on Korean coasts The TCS haplotype network analysis of 229 A. Table 2. Genetic diversity indices for the COI of 229 polyp samples of Aurelia coerulea collected from four regions (R1, R2, R3, and R4) on western and southern coasts of Korea Region N Nh Ps D FShπ R1 79 21 27 -0.81 -5.215 0.8079 ± 0.0265 0.005316 ± 0.002973 R2 60 20 27 -1.284 -7.152* 0.8554 ± 0.0349 0.004577 ± 0.002626 R3 36 15 22 -0.654 -3.301 0.9 ± 0.0348 0.005647 ± 0.003183 R4 54 21 27 -1.05 -7.541* 0.8819 ± 0.032 0.00529 ± 0.002979 Total 229 53 51 -1.487 -25.442* 0.8788 ± 0.0169 0.005532 ± 0.003052 N, Number of sample; Nh, Number of haplotype; Ps, Number of polymorphic site; D, Tajima’D score; FS, Fu’s fs score; h, Haplotype diversity; π, Nucleotide diversity; *, p < 0.05. Table 3. Analysis of molecular variance (AMOVA) results showing partitioning of genetic variation of (A) all populations and (B) between two haplogroups (A and B) Source of variation d.f. Sum of Squares Variance Components % of Variation Fixation indices P-value (A) All regions Among regions 3 37.792 -0.06 -2.78% FCT = -0.028 Not significant Among locations Within regions 5 46.376 0.409 18.94% FSC = 0.184 p = 0 Within Populations 220 398.667 1.812 83.84% FST = 0.162 p = 0 (B) Haplogroups (A and B) Among haplogroups 1 278.06 2.673 77.35% FST = 0.774 p = 0 Within Populations 216 169.066 0.783 22.65% page 5 of 13Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan Fig. 2. Pairwise FST data based on the 53 COI haplotypes of Aurelia coerulea analyzed between nine locations, within four regions (R1, R2, R3 and R4) of Korean coasts. The abbreviations of the locations are listed in table 1. Fig. 3. Haplotype (TCS) network showing the phylogeography of the 53 COI haplotypes of Aurelia coerulea collected from western and southern coasts of Korea. Each circles represent individual haplotypes and the size of the cycle is the proportion of haplotype frequency. The color of each segment indicates the geographic origin of that fraction of the samples. Each line between haplotypes represents a 1-nucleotide mutational change. Small black dots indicate unsampled haplotypes. JG DJ DS BJ DC GL UY ND GS JG DJ DS BJ DC GL UY ND GS Pairwise FST COI R1 R2 R3 R4 0.1 0 0.1 0.2 0.3 0.4 0.5 FST 1 71 23 10 - 13 2 - 4 6 - 8 R1 R2 R3 R4 H09 H06 H20 H04 H14 H10 H03 H26 H21 H19 H12 H30 H49 H47 H48 H22 H18 H07 H29 H23 H08 H40 H37 H43 H25 H28 H52 H50 H41 H02 H45 H34 H38 H39 H32 H33 H17 H15 H05 H13 H16 H35 H27 H31 H24 H36 H42 H46 H44 H01 H51 H53 H11 Haplogroup A Haplogroup B page 6 of 13 Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan coerulea polyps revealed the phylogenetic similarities and geographical distribution of COI haplotypes in Korean coastal waters. Overall, the 53 haplotypes were divided into two distinct haplogroups (A and B) (Fig. 3). Haplogroup A was formed by 19 haplotypes; all the haplotypes in the group were linked to the most frequent haplotype (H04), showing a star-like phylogeny consistent with a pattern of population expansion. In contrast, haplogroup B comprised 24 haplotypes that formed a branched haplotype network. The two haplogroups were separated by five mutational steps, and there were 10 intermediate haplotypes that did not belong to either group (H02, H15, H17, H32, H33, H34, H38, H39, H41, and H45). Overall, the most frequent haplotype was H04 from haplogroup A (71 individuals), and this haplotype was shared among all four regions. Within haplogroup B, H09 was the most frequent haplotype. In total, there were 13 haplotypes that shared more than one region (Fig. 3), and the rest were unique to a single region. H26 and H06 were the most frequent unique haplotypes found in haplogroup A and B, respectively. AMOVA analysis based on the two haplogroups (A and B) yielded a significant genetic variation, with 77.35% of the genetic differentiation occurring among populations and 22.65% occurring within populations. The FST value was 0.774, highlighting a gene flow restriction between the two haplogroups (Table 3). The 53 haplotypes were distributed across the four regions of Korean coastal waters. The haplotypes belonging to the two haplogroups coexisted in all four regions, despite significant genetic variation between the two haplogroups (Fig. 4). As mentioned above, haplotype H04 was the most frequently found haplotype in all four regions. However, in R1, H09 and H06 were codominant with H04 (Fig. 4). Analysis of the PCA was carried out to summarize the genetic diversity on a geographical basis. The analysis revealed the distribution of all 53 haplotypes in the four regions (Fig. 5). Each green arrow represents the geographical locations, with angles between arrows indicating the significance of the correlation and the distance between the dots in the diagram approximating the dissimilarity of the haplotypes in the four regions. Most of the haplotypes (Black dots) are less frequent and have low variation in the regions. Haplotypes H09 and H06 in R1 and H26 in R2 showed the greatest variation, while H04 is the most frequent. Fig. 4. Distribution of the COI haplotypes of Aurelia coerulea in four regions (R1, R2, R3 and R4) in the Korean coastal waters. The small circles represent the proportion of haplogroup A, haplogroup B and the intermediate haplotypes in each region. The big circles represents the proportion of the haplotypes found in each region, and the number in the middle indicates the total number of individuals sampled in each region. Each color in the circles represents the frequencies of the dominant haplotype. Others are the total frequencies of the haplotypes that are not shown in the respective circles. H04 H06 H09 79 H04 H10 H11 H14 H26 60 H04 H03 H10 H20 H36 36 H03 H04 H20 54 R1 R2 R3 R4 Haplogroup A Haplogroup B Intermediate haplotypes page 7 of 13 Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan Phylogenetic analysis The present phylogenetic tree was constructed based on the neighbor-joining (NJ) method using the COI gene sequences of our 53 haplotypes and 115 COI sequences of the global A. coerulea population retrieved from the GenBank database (ESM 1, Table S2 and S3). The tree showed that all global A. coerulea populations were divided into two main clades that were similar to the TCS haplotype network of the Korean populations. All haplotypes in haplogroup A from the TCS network clustered in Clade I (Fig. S1A), while all the haplotypes in haplogroup B clutered in Clade II (Fig. S1B). Haplotypes from Japan, China, Europe, and Australia nested in both clades. However, most of the European and Australian haplotypes nested in Clade I. In addition, American haplotypes clustered only in Clade I. The intermediate haplotypes in TCS network were indistinguishable and clustered in Clade I of the phylogenetic tree. DISCUSSION Genetic structure In the present study, we found high genetic diversity in all regions, especially in R3 and R4. Overall, high genetic diversity was found in Korean coastal waters, which was similar to a previous finding in Japanese populations (h = 0.87 ± 0.06; π = 0.0063 ± 0.0037) (Dawson et al. 2005). The high genetic diversity found in this study was expected because Aurelia coerulea is believed to be a native species in the western North Pacific (Dawson et al. 2005). In addition, Fofonoff et al. (2018) mentioned that the species is most likely a native from southern Korea to northern Japan and possibly China. Usually, high genetic diversity is observed in the native populations and is favoured by a large stable population size and positive demographic balance (Grant and Bowen 1998; Micheli et al. 2015), which improves the resistance and resilience of populations to environmental changes (Roger et al. 2012). The genetic data presented here support the idea that A. coerulea is likely native to the western North Pacific, as suggested by Dawson et al. (2005). However, the negative and significant values shown by Tajima’s D and Fu’s FS tests indicated that there has been a recent population expansion and an increase in effective population size in Korean coastal regions. The AMOVA analysis based on all A. coerulea populations showed that most of the genetic differentiation occurred within populations. This implies that there was no restricted gene flow between the four regions (except between R1 and R2) and shows a sympatric distribution of the two haplogroups. Restricted gene flow is caused primarily by geographical barriers and low dispersal ability of species (Liu et al. 2019). Pairwise FST analysis showed that most of the genetic differentiation occurred between region R1, specifically DS, and the rest of the regions. This is apparently due to haplotypes H06 and H09, which were mainly found in DS (Fig. S2). Similarly, haplotype H26 was restricted to R2, resulting in unequal gene flow and high genetic differentiation between R1 and R2 despite a short geographic distance. This is clearly shown in the PCA analysis. Usually, genetic differentiation between more distant populations is higher due to substantial differences in allele frequencies (Liu et al. 2019), which could be interpreted as isolation-by-distance or allopatric speciation (Coyne and Orr 2004; Zhan et al. 2009). However, in this case, it could not be isolationby-distance; rather, R1 and R2 displayed a divergent pattern similar to parapatric speciation (Butlin et al. 2008). That is, there was restricted gene flow between the diverging populations. The reason for this unique differentiation between R1 and R2 may be due to hydrodynamic circulation, habitat fragmentation formed by marine gyres and currents, or even anthropogenic activities (Otto et al. 1990; Launey et al. 2002; Zhan et al. 2009; Dong et al. 2015; Liu et al. 2019). Fig. 5. PCA plot showing the distribution of Aurelia coerulea haplotypes in the four regions (R1, R2, R3, and R4) of Korean coastal waters. The black dots are the haplotypes that are less frequent and have less variation in the four region. The coloured dots are haplotypes that are more frequent and have more variation between the four regions. The green arrows represent parameters pointing in the direction of the four geographical locations (R1, R2, R3, and R4). R1 R2 R3 R4 -10 -5 0 5 10 15 20 25 30 35 PC1 (78.593%) -20 -15 -10 -5 0 5 10 PC2 (15.446%) H04 H09 H06 H26 Other haplotypes page 8 of 13Zoological Studies 60:63 (2021) © 2021 Academia Sinica, Taiwan Phylogeography and intra-specific variation in A. coerulea Our TCS haplotype network analysis revealed a well-resolved phylogeography of A. coerulea in Korean coastal waters. The 229 polyp samples collected were clearly divided into two major haplogroups (A and B) separated by five mutational steps with intermediate haplotypes. AMOVA constructed based on these two distinct haplogroups showed high genetic differentiation among populations with a high FST value of 0.77353. This clearly indicates high genetic variation and restricted gene flow between the two groups. The haplogroups displayed different patterns, with haplogroup A having a star-like haplotype network, providing support for a recent population expansion, as described by Avise (2000) and Froufe et al. (2016). Haplogroup B showed a branched haplotype network, which is thought to be the nature of a stable and older clade (Ferreri et al. 2011). Interestingly, the haplotypes of the two haplogroups coexisted in all regions in Korean coastal waters, and the coexistence of these divergent haplogroups in the same geographical locations suggests sympatric development (Avise 2000). Population subdivision and intraspecific phylogeographic structures are commonly observed in jellyfish (e.g., Holland et al. 2004; Dawson 2005; Ramšak et al. 2012; Lee et al. 2013) and other organisms, including bluegill sunfish (Lepomis macrochirus; Avise et al. 1984), rock slaters (Ligia occidentalis; Markow and Pfeiler 2010), and freshwater mussels (Potomida, Unionida; Froufe et al. 2016). However, these are all examples of allopatric differentiation, the most commonly encountered situation. It has been rare to observe large mtDNA differences within the same species collected at any given geographic site (Avise et al. 1987). Though this large mtDNA difference was recently reported in the jellyfish Pelagia noctiluca (Stopar et al. 2010; Glynn et al. 2016), we we did not find anything like this reported in A. coerulea. Yi and Kim (2020) found similar significant genetic differentiation among ascidian Ciona savignyi populations along Korean coasts. They suggested that a possible northward post-glacial colonisation scenario in the western North Pacific around 15 thousand years ago left evidence of the expansion and differentiated genetic structure of C. savignyi populations and other coastal marine organisms in Korea (Yi and Kim 2020). We did not estimate molecular clock calibration for scyphozoans in the present study, because a reliable molecular clock calibration for scyphozoans is not available (Holland et al. 2004); therefore, the estimated separation time could not be determined. However, with high genetic variation and restricted gene flow found between the two haplogroups in this study, the above hypothesis gives a possible explanation for the COI variation found in A. coerulea. Based on this and inferences made by several previous studies, Korean coastal regions may be considered a zone of secondary admixture between allopatrically evolved populations (Avise et al. 1987; Avise et al. 1998; Avise 2000; Hewitt 2000; Beebee and Rowe 2004; Provan and Bennett 2008). In addition, for the secondary contact and expansion event, haplogroup B, with its diverse network topology, might be the older clade that settled a long time ago in Korean coastal waters. However, the starlike phylogeny of haplogroup A is generally interpreted as indicative of a population that recently expanded (Slatkin and Hudson 1991; Avise 2000; Froufe et al. 2016). This is supported by the overall neutrality test from this study. Therefore, there is strong evidence for past population expansion through introductions by anthropogenic means or natural dispersal, shaping the genetic structure of A. coerulera in Korean coastal waters. Global distribution of A. coerulea As mentioned above, the moon jellyfish A. coerulea has the broadest geographic range of all species in the genus Aurelia to date (Ki et al. 2008; Scorrano et al. 2017). The species was suggested to be native to the western North Pacific (Japan, Korea, and China) and dispersed globally from Japan by anthropogenic means (Dawson et al. 2005). This is supported by the present phylogenetic analysis based on the global A. coerulea population, which showed that haplotypes in haplogroup A and B are widely distributed globally, and the most frequent haplotype was clustered with Chinese and Japanese haplotypes. Due to the overlapping distribution of the two haplogroups in Korean coastal regions, it is difficult to determine the geographical origin of the two haplogroups. However, the present data highlight the possibility that, after the speciation events, haplogroup B became largely concentrated in the Bohai Sea and the Yellow Sea, while haplogroup A may have concentrated in Japanese coasts and dispersed naturally, mixing locally in the western North Pacific. However, dispersal beyond this region to other global locations occurs through anthropogenic means (Dawson et al. 2005). The presence of shared haplotypes among distant locations highlights recent connections through human-mediated dispersal (Lejeusne et al. 2014). page 9 of 13Zoological Studies 60:63 (2021)