Insights Into the Environmental Impact on Genetic Structure and Larval Dispersal of Crown-of-Thorns Starfish in the South China Sea

The coral-eating crown-of-thorns starfish (COTS; Acanthaster spp.) play a major role in coral reef degradation in the Indo-Pacific region. However, the impacts of environmental factors on the phylogenetic and genetic characteristics of COTS in the northern Indo-Pacific convergence region remains unclear. We used mitochondrial DNA (mtDNA) and microsatellite markers to analyze the phylogenetic relationship, demographic history, genetic diversity and genetic structure of COTS in the South China Sea (SCS) and explored the impact of environmental factors on historical population expansion, genetic differentiation and larval dispersal. There was a clear signature of a population expansion in the SCS using the mtDNA marker. According to microsatellite loci analysis, COTS have high genetic diversity in the SCS. STRUCTURE analysis indicated that COTS in the Pacific Ocean can be divided into four subgroups: the SCS, Western Pacific, Pacific equatorial current affected zone, and Pacific insular atolls populations in the Pacific Ocean. Fst-statistical analysis revealed positive correlations between the Fst values and geographic isolation for all sampling sites. Additionally, there were no clear associations between the Fst values and chlorophyll a concentrations among coral reefs in the SCS; however, there were significant positive associations between the Fst values and particulate organic carbon (POC) concentrations within small geographic distances. These results suggest that COTS underwent historical population expansion after the Last Glacial Maximum, possibly followed by coral population expansion. The genetic structure of COTS populations may have been shaped by distinct nutrient concentrations, particularly those of POC, over small geographic distances. Moreover, ocean currents provide a potential dispersal mechanism for COTS larvae in the SCS. This study demonstrates that environmental and oceanographic factors play important roles in shaping the genetic characteristics and larval dispersal of COTS populations in the northern Indo-Pacific convergence region.

The coral-eating crown-of-thorns starfish (COTS; Acanthaster spp.) play a major role in coral reef degradation in the Indo-Pacific region. However, the impacts of environmental factors on the phylogenetic and genetic characteristics of COTS in the northern Indo-Pacific convergence region remains unclear. We used mitochondrial DNA (mtDNA) and microsatellite markers to analyze the phylogenetic relationship, demographic history, genetic diversity and genetic structure of COTS in the South China Sea (SCS) and explored the impact of environmental factors on historical population expansion, genetic differentiation and larval dispersal. There was a clear signature of a population expansion in the SCS using the mtDNA marker. According to microsatellite loci analysis, COTS have high genetic diversity in the SCS. STRUCTURE analysis indicated that COTS in the Pacific Ocean can be divided into four subgroups: the SCS, Western Pacific, Pacific equatorial current affected zone, and Pacific insular atolls populations in the Pacific Ocean. Fst-statistical analysis revealed positive correlations between the Fst values and geographic isolation for all sampling sites. Additionally, there were no clear associations between the Fst values and chlorophyll a concentrations among coral reefs in the SCS; however, there were significant positive associations between the Fst values and particulate organic carbon (POC) concentrations within small geographic distances. These results suggest that COTS underwent historical population expansion after the Last Glacial Maximum, possibly followed by coral population expansion. The genetic structure of COTS populations may have been shaped by distinct nutrient concentrations, particularly those of POC, over small geographic distances. Moreover, ocean currents provide a potential dispersal mechanism for COTS larvae in the SCS. This study demonstrates that environmental and oceanographic factors play important roles in shaping the genetic characteristics and larval dispersal of COTS populations in the northern Indo-Pacific convergence region.

INTRODUCTION
Coral reefs are one of the most biologically diverse ecosystems in the oceans (Reaka-Kudla, 1997;Robert et al., 2002;Wilkinson, 2008). However, coral reef ecosystems have been severely threatened by anthropogenic global climate change and this has become especially apparent since the 1970s (Silverstein et al., 2015). Climate change not only leads to drastic declines in coral cover (Hughes et al., 2018;Stuart-Smith et al., 2018), but also appears to promote outbreaks of coral predator crown-of-thorns starfish (COTS; Acanthaster spp; Uthicke et al., 2015a). COTS are the most important biological threat to coral reef ecosystems in the Indo-Pacific region (Pratchett et al., 2017a), and COTS outbreaks not only lead to substantial loss of live coral cover, but also impair the integrity and resilience of reef ecosystems (Death et al., 2012;Kayal et al., 2012;Mellin et al., 2019). In the Indo-Pacific region, at least 246 outbreaks of COTS have been recorded since 1990, which is triple the number of population explosion events reported prior to 1990 . Additionally, coral reefs in the Pacific region experience a high frequency of COTS outbreaks, including the Great Barrier Reef (GBR; Death et al., 2012;Pratchett et al., 2014), and those in Japan (Yamaguchi, 1986;Suzuki et al., 2012;Yasuda, 2018), French Polynesia (Kayal et al., 2012), Guam (Chesher, 1969), and other regions, which has resulted in high levels of coral mortality and reduced live coral cover. Although many COTS outbreak events have been observed, the causes of dramatic COTS population explosions remain equivocal. It is still unclear whether outbreaks start from a single reef or arise simultaneously from separate locations due to limited temporal and spatial resolution of monitoring (Pratchett, 2005).
Crown-of-thorns starfish have high fecundity, and a female COTS can produce 30-70 million eggs (resulting in about 10 million fertilized eggs) per spawning period (Birkeland and Lucas, 1990;Babcock et al., 2016;Caballes et al., 2021;Pratchett et al., 2021). Thus, many studies have suggested that COTS outbreaks can be attributed to higher survival rate during early life stages, and that a small increase in the survival rate of the COTS larvae could lead to rapid increases in population size (Birkeland and Lucas, 1990). Elevated food supply (Brodie et al., 2005), anthropogenic eutrophication of seawater (Lucas, 1973;Birkeland, 1982;Fabricius et al., 2010), reduced predation pressure due to overfishing (Cowan et al., 2016) have been associated with increased larval survival. Notably, the geographical dispersal of COTS plays an important role in the formation of secondary outbreaks of COTS (Pratchett, 2005;Yasuda et al., 2009;Wooldridge and Brodie, 2015;Tusso et al., 2016), given that the duration of the planktonic larval stage (PLS) ranges between 9 and 42 days . Weak or variable longshore currents may promote strong larval retention or limited dispersal, which are likely to give rise to primary outbreaks (Pratchett, 2005;Wooldridge and Brodie, 2015). However, strong directional longshore currents increase the likelihood of inter-reef dispersal, which leads to outbreaks once they become established (Wooldridge and Brodie, 2015). The multiple factors that lead to increased COTS larval survival rates remain unclear, and they may show spatial differences among distinct areas in the Indo-Pacific region. In addition, the dispersal routes and spatial distribution of COTS larvae may also vary between sites (Uthicke et al., 2015b;Yasuda et al., 2015a;Suzuki et al., 2016). Therefore, the knowledge of relationships between environmental-oceanographic factors and population size, demographics and structure of COTS population need to be improved.
The survivorship, settlement, and development of COTS larvae has been difficult to document in the wild, so direct observations have proven ineffective for COTS outbreaks, population demographics and dispersal events at different life stages. Model-based statistics of population genetics have increasingly been used as an indirect and effective method to infer the patterns of historical demographics, genetic diversity, genetic structure and migration routes for COTS, which could aid in understanding the relationship between molecular ecological characteristics of COTS and environmentaloceanographic factors among distinct coral habitats. Vogler et al. (2008) identified four clades of COTS were within the Indo-Pacific regions based on mitochondrial DNA (mtDNA) gene markers, which included Pacific, Red Sea, Southern Indian Ocean and Northern Indian Ocean. These clades were estimated to have diverged 1.95-3.56 Mya in the Pliocene-Early Pleistocene period. Unexpectedly, results from mtDNA studies have shown that genetic differentiation of COTS populations have clear relationships with regions. For instance, Timmers et al. (2012) showed that there are shared mitochondrial haplotypes between the South Central and Northwest Pacific and that their haplotypes do not strictly cluster according to geographic region. Vogler et al. (2013) also found a large geographic cluster in Western Pacific localities, with shared haplotypes across the entire range from the GBR to the Philippines. In addition, studies based on sensitive microsatellite loci makers revealed that COTS populations are genetically differentiated, there is contemporary dispersal among distinct coral reefs, and there are established relationships between ocean currents and maximum larval dispersal distance (Yasuda et al., 2009Tusso et al., 2016;Harrison et al., 2017). A large genetic break has been documented between the Indian Ocean and the Pacific Ocean corresponding to mtDNA study by Vogler et al. (2008), and there is significant genetic differentiation of COTS populations among Southeast Africa, Northwestern Pacific, Palau, the North Central Pacific, GBR, Fiji, and French Polynesia (Yasuda et al., 2009. Furthermore, it has also been suggested that western boundary currents played an important role in the dispersal of COTS larvae in the Southeast African, Northwestern Pacific, and GBR oceanic regions (Yasuda et al., 2009). Tusso et al. (2016) also found that COTS populations show significant differentiation among Pacific insular atolls (e.g., Moorea and Johnston Atoll), and there was an effective dispersal pattern and gene flow of the COTS population in the Pacific equatorial current affected zone (e.g., Guam, Kingman reef and Swains; Tusso et al., 2016). However, in previous population genetics studies of COTS, generally conducted in areas within or near the open ocean (Yasuda et al., 2009Vogler et al., 2012Vogler et al., , 2013Tusso et al., 2016), strong ocean currents appear to facilitate effective larval dispersal of marine organisms. There are very few studies about the associations between genetic differentiation and connectivity of COTS populations in northern Indo-Pacific convergence region. Importantly, this region has special geographic barriers (e.g., mainland and archipelago) that may form genetic barriers against the mixing of COTS populations between marginal seas and the open ocean.
The South China Sea (SCS) spans the tropical and subtropical climate zones in the northern Indo-Pacific convergence region, and contains numerous coral reefs and communities across 19 degrees of latitude (Spalding, 2001;Yu, 2012;Chen et al., 2020). The abundant growth of corals in SCS provides favorable habitat for COTS populations (Li et al., 2019). Moreover, the SCS is connected to the Indian Ocean and the Pacific Ocean via several straits, with seawater exchange through the Luzon Strait, Taiwan, and Mindoro (Yu et al., 2010). Some studies have shown that there has been significant genetic differentiation and connectivity between reef organisms in the SCS and the Pacific Ocean (e.g., Symbiodiniaceae), which may be due to large biogeographic barriers and latitudinal environmental factors (Chen et al., 2020). In addition, reef around the Paracel Islands (Xisha Islands) had experienced dramatic population explosions of COTS that destroyed > 95% of live coral cover from 2005 to 2010 (Li et al., 2019). Unexpectedly, in less than 8 years, COTS outbreaks occurred in the Paracel Islands in 2019 again (Li et al., 2019). Nevertheless, the knowledge about genetic structure, connectivity, larval dispersal and its environmental impact factors of COTS in the SCS is still scant.
In this study, mtDNA and microsatellite genetic markers were used to analyze the phylogenetic relationship, demographic history, genetic diversity, genetic structure and differentiation, and recent gene flow of the COTS populations in the SCS. Previously published population genetic datasets of COTS in the Indo-Pacific region were also used for analyses and to draw comparisons (Vogler et al., 2008;Tusso et al., 2016). These analyses improve our understanding of COTS historical population expansion, long-distance larval dispersal, and genetic differentiation of COTS between SCS and other coral habitats.

Samples Collection and DNA Extraction
In total, 163 adult COTS with similar sizes (ranging from 210 to 220 mm) were collected by SCUBA diving in June 2018, which ensure that individual samples were approximately the same age (2-3 years old; Pratchett et al., 2014). The sampling sites included five distinct coral reefs in the SCS: First Thomas Shoal (FTS; Xinyi Reef), Passu Keah (PK; Panshiyu), Drummond Island (DI; Jinqing Island), Bombay Reef (BR; Langhua Reef), and Pattle Island (PI; Shanhu Island) (Figure 1 and Table 1). Tube feet of COTS samples were cut off using sterile medical operation scissors. All tube feet were transferred into 50-mL cryogenic vials (Corning, Mexico), 95% ethanol was added, and tube feet samples were stored at -20 • C until DNA extraction. DNeasy R Blood and Tissue Kits (Qiagen, Hilden, Germany) were used to extract total DNA from tube feet samples, following the manufacturer's instructions.

Phylogenetic Tree Construction and Demographic History
A fragment of the mitochondrial cytochrome oxidase subunit I gene (COI) was amplified using COTS_COI_R4734 5 -GCCTGAGCAGGAATGGTTGGAAC-3 and COTS_COI_ R5433 5 -CGTGGGATATCATTCCAAATCCTGG-3 (Vogler et al., 2008). PCR was conducted with 25 µL Premix Taq TM (TaKaRa Taq TM Version 2.0 Plus dye, TaKaRa, Dalian, China), ∼75 ng of DNA, 2 µL primer, and ddH 2 O to a total volume of 50 µL. Based on ABI GeneAmp 9700 thermal cycler, PCR reactions were carried out under the following conditions: 2 min at 94 • C, followed by 35 cycles of 94 • C for 20 s, 59 • C for 45 s, 72 • C for 1 min, and a final extension at 72 • C for 10 min. Amplification products were directly sequenced on an ABI 3730XL DNA Analyser (Sangon Biotech Company, Shanghai, China). COI sequences were then aligned using the ClustalW algorithm implemented in MEGA 6. Phylogenetic trees of COTS sequences were inferred using maximum likelihood (ML) and neighbor-joining (NJ) methods. ML analyses for sequence data were performed with a K2 + G model. The robustness of the ML and NJ trees was tested with 1000 bootstrap replicates. Novel COI sequences of COTS in the SCS can be accessed in GenBank under accession numbers MW238184-MW238343. In addition, DnaSp6 was used to analyze the number of haplotypes (H), haplotype diversity (h), the number of effective haplotypes [EH; HE = 1/(1 − h)], and nucleotide diversity (π), and NETWORK v10.1.0.0 was used to visualize haplotype variability and clustering based on the median-joining method. To understand the demographic history of COTS in the SCS and investigate signatures of historical population expansions, Fu's Fs (Fu, 1997) and Tajimas's D (Tajima, 1989) were calculated using Arlequin v 3.5 (45,000 replicates; Excoffier and Lischer, 2010).

Microsatellite Analysis
Nine polymorphic microsatellite loci were used to analysis COTS genotypes (Supplementary Table 1; Yasuda et al., 2006Yasuda et al., , 2007Tusso et al., 2016). Each locus was amplified in 50 µL reaction volumes containing ∼70 ng DNA, 2 µm primer, 24 µL of TaKaRa Premix Taq TM and ddH 2 O. The PCR were performed under the following conditions: 94 • C 2 min, followed by 35 cycles of 94 • C for 20 s, 51 ) for 20 s, and a final extension at 72 • C for 5 min. ABI 3730XL Genetic Analyzer was used to analyze PCR products at the Sangon Biotech Company (Shanghai, China). Sizing was achieved using LIZ500 internal standards. The GeneMapper v.3.2 was used to identify allele sizes. To explore the genetic differentiation and gene flow of COTS populations between the SCS and the Pacific Ocean, we combined microsatellite data of COTS from the Pacific Ocean (Western Pacific, Pacific equatorial current affected zone, and Pacific insular atolls), which were published by Tusso et al. (2016). The microsatellite genotype data sets have been deposited in Supplementary Dataset 1.
Hardy-Weinberg equilibrium (HWE) and linkage disequilibrium exact tests were conducted by Genepop v4.3.1. Myriad v1.2 software was used to test for the presence of null alleles using Bonferroni-adjusted 95% confidence intervals (Antonio, 2018). The descriptive statistics included the number of alleles (Na), observed heterozygosity (Ho), expected heterozygosity (He), Private alleles (Pa) and Information index (I), which were conducted by GeneAlex v.6.5 (Peakall and Smouse, 2010). Inbreeding coefficient (Fis) was calculated by Genepop v4.3.1, and allelic richness (Ar; with 2000 randomizations) and was used to determine differences in genetic diversity using FSTAT v.2.9.4 (Goudet, 2001). A General Linear Model (GLM) analysis was used to assess differences in genetic diversity with 1000 bootstrap permutations, and the groups of ocean regions was regarded as a fixed factor. The Student-Newman-Keuls (SNK) test was used for post hoc multiple comparisons of significant GLM test results. STRUCTURE v2.4 was used to analyze the genetic clusters, and Markov chain runs consisted of an initial "burn-in" of 250,000 steps followed by a final 1,000,000 iterations (Evanno et al., 2010). Fifteen independent runs were performed for 1-19 Ks. The most likely value for K was tested by plotting the log probability of data over multiple runs and comparing that with delta K as implemented in Structure Harvester v.0.6.94 (Earl and VonHoldt, 2012), and the results of Ln Pr(X\K) and delta K are shown in the supporting information (Supplementary Figure 1). To choose the appropriate K value the results were averaged using CLUMPP v.1.1.2 (Jakobsson and Rosenberg, 2007), and genetic structure plots were manipulated using DISTRUCT (Rosenberg, 2010). Analyses of molecular variance (AMOVA) were conducted on geographic populations with regional structuring as defined by STRUCTURE and by sampling sites to determine the degree of differentiation, all of which was performed in GeneAlex v.6.5 (Peakall and Smouse, 2010). Genetix software was employed to calculated pairwise Fst values, using Weir and Cockerham genetic distance with 10,000 permutations to determine significance. Correlations between  (Bohonak, 2002).

Environmental Characteristic Analysis
The genetic differentiation and long-distance dispersal of COTS larvae may be closely related with geographic isolation (Dight et al., 1990;Vogler et al., 2012), nutrient concentrations (Fabricius et al., 2010;Pratchett et al., 2017a,b) and ocean currents (Yasuda et al., 2009;Vogler et al., 2012;Tusso et al., 2016). As the food source of COTS larvae was diverse and complicated (Ayukai, 1994;Okaji et al., 1997;Nakajima et al., 2016;Mellin et al., 2017), the chlorophyll a (Chl a, µg·L −1 ) and particulate organic carbon (POC, µg·L −1 ) concentrations were collected together as proxies for the nutrient concentration, which can represent the phytoplankton biomass and particulate organic matter concentrations in the sampled coral reefs and surrounding sea water, respectively. Aqua-MODIS monthly average Chl a and POC datasets from 2016 to 2020 were obtained from the NASA Giovanni satellite 1 . Additionally, correlations between Chl a and Fst and between POC and Fst were tested by Mantel tests using IBD software v.1.5.2 (Bohonak, 2002)

Phylogenetic Analyses and Demographic History of Crown-of-Thorns Starfish Populations
The phylogenetic analysis of COTS mtDNA showed that almost all COI sequences of COTS in the SCS do not show obvious phylogenetic divergence with Pacific populations (Figure 2A).
In addition, COTS populations in SCS also have significant phylogenetic divergence with those in the Red Sea, Southern Indian Ocean and Northern Indian Ocean (Figure 2A). For COTS populations within the SCS, there were no significant distinct genetic clusters in the median-joining haplotype ( Figure 2B; separated by at least 18 mutation steps; Vogler et al., 2013). This result also showed no clear association between haplotypes and geographic factors, and two dominant haplotypes were widely shared among the populations of the coral reefs in the SCS ( Figure 3B). The haplotype diversity in COTS populations in the SCS was low, and only the h of COTS population in BR was >0.1 (h = 0.910; Table 1). Additionally, low nucleotide diversity was determined, and π of all populations was <0.01 in the SCS (Table 1). However, there was a clear signature of a population expansion identified by mismatch distribution analysis ( Figure 2C) and the significance of Fu's Fs and Tajima's D ( Table 1).

Hardy-Weinberg Equilibrium and Linkage Disequilibrium
For microsatellite loci analysis, based on the multi-locus test for five COTS populations in the SCS, most microsatellite loci were at or near HWE after the Bonferroni correction test (p < 0.05). In addition, there were some deviations in COTS populations from the SCS; for example, Aya02 and Yukina06 deviated from HWE in the FTS, PK, and BR. However, almost all of these loci exhibited HWE in the COTS populations of the SCS. Moreover, none of the pairwise comparisons between loci within COTS populations from the SCS and the Pacific Ocean exhibited significant evidence of linkage disequilibrium after Bonferroni correction. Therefore, these loci were included in the data set for subsequent analyses.

Genetic Diversity
The results of the genetic diversity analysis using the microsatellite gene marker were not consistent with those of the mtDNA sequences analysis. COTS populations in the SCS, Western Pacific, and Pacific equatorial current affected zone had higher average allelic richness than that in the Pacific insular atolls (GLM; F = 11.260, p < 0.01; Figure 3A). The average allelic richness was highest in the SCS (Ar mean = 3.829; SD = 0.763), followed by the Pacific equatorial current affected zone (Ar mean = 3.603; SD = 0.716), Western Pacific (Ar mean = 3.538; SD = 0.675), and Pacific insular atolls (Ar mean = 2.478; SD = 1.002). The COTS population in DI had the highest average allelic richness (Ar mean = 3.914; SD = 0.803), and the lowest average allelic richness was in the Johnston Atoll COTs population (Ar mean = 2.220; SD = 1.153). In addition, the observed heterozygosity of COTS populations in the SCS and Pacific equatorial current affected zone were higher than those in Western Pacific and Pacific insular atolls (GLM; F = 31.441, p < 0.01; Figure 3B). However, the changing rule of the expected heterozygosity was consistent with average allelic richness among distinct regions, and the COTS populations in the SCS, Western Pacific, and Pacific equatorial current affected zone regions had higher expected heterozygosity than that in the Pacific insular atolls (GLM; F = 19.664, p < 0.01; Figure 3C). It is worth noting that the COTS populations in the SCS have relatively lower inbreeding coefficients (Fis) and higher private alleles (Pa) than that in other regions (Table 1). Accordingly, COTS populations in the SCS, Western Pacific, and Pacific equatorial current affected zone have higher genetic diversity than that in the Pacific insular atolls.

Genetic Differentiation and Structure
The analysis by STRUCTURE identified two peaks (K = 2 and K = 6) as the most likely number of COTS genetic clusters based on estimates of delta K and LnP(K) (Figures 4A,B and Supplementary Figure 1). At K = 2, the Pacific equatorial current affected zone formed a separate cluster from the other COTS samples (Figure 4A). At K = 6, COTS populations from the SCS, Western Pacific, Pacific equatorial current affected zone, and Pacific insular atolls comprised separate cohesive genetic units ( Figure 4B). COTs populations in the Pacific equatorial current affected zone had admixed genotypes, and two main clusters (orange and dark blue) were identified in samples from Guam, Kingman Reef and Swains. The samples from Japan, the Philippines, the GBR and Vanuatu belonged mostly to the same inferred population (green), which was located in the Western Pacific. An additional genetic cluster (bright yellow) was found at low frequencies but with high individual membership probabilities in the Johnston Atoll ( Figure 4B). Interestingly, COTS samples from the SCS formed two separate genetic clusters, namely the First Thomas Shoal-Passu Keah cluster (FTS-PK; blue) and the Drummond-Bombay Reef cluster (DI-BR; gray; Figure 4B). In addition, COTS populations of PI and Moorea had complex compositions of genetic clusters, which may be associated with insufficient sample sizes. Thus, microsatellite loci data sets from these two locations were not used to calculate Fst values.
Fst values were significant in 75% of all pairwise comparisons (with 10,000 permutations), with larger values obtained between distinct regions and between more distinct sites ( Table 2). There were mainly insignificant genetic differentiation rates among COTS populations within the same region. For instance, FTS and PK in the SCS (Fst = 0.001, FDR, p > 0.05), the Philippines and the GBR in the Western Pacific (Fst = 0.041, FDR, p > 0.05), and some Guam sampling sites and Kingman Reef in the Pacific equatorial current affected zone (Fst values from 0.008 to 0.039, FDR, p > 0.05). AMOVA also found that there were significant genetic differentiation rates among regions ( RT = 0.129, p < 0.001, Table 3), among sampling sites ( PR = 0.066, p < 0.001, Table 3) and within COTS individuals ( PT = 0.187, p < 0.001, Table 3). It was noted that the SCS, Western Pacific, Pacific equatorial current affected zone and Pacific insular atolls yielded highly significant values (percent of variation: 13%), which indicated strong geographical differentiation among the four regions. IBD results also indicated significant positive correlations between the Fst values and geographic distances (Mantel test, Randomization = 10,000, R 2 = 0.1347, p = 0.0099; Figure 5A), which was consistent with the results of the AMOVA. It is worth noting that nutrient concentration was significantly distinct among coral habitats within short geographic distances inside the SCS (GLM; Chl a: F = 92.715, p < 0.001; POC: F = 140.043, p < 0.001; Figure 6). Additionally, the results of the Mantel test and DistLM found that there was no significant correlation between the Fst values and Chl a (Mantel test, Randomization = 10,000, R 2 = 0.3121, p = 0.1150; DistLM, pseudo-F = 1.815, p = 0.2492; Figure 5B); however, the Fst values were significantly related to POC (Mantel test, Randomization = 10,000, R 2 = 0.8670, p = 0.0010; DistLM, pseudo-F = 26.07, p = 0.0070; Figure 5C) within the SCS. Thus, nutrient concentration, especially POC concentration, affect the genetic differentiation of COTS populations among distinct coral habitats, which only had slight geographic isolation from each other.

Historical Population Expansion of Crown-of-Thorns Starfish May Be Associated With Coral Range Expansion After the Last Glacial Maximum in the Pleistocene
The mtDNA marker (4734-5433) for phylogenetic analysis of COTS populations did not belong to the highly variable region   of COI fragment (Yasuda et al., 2009;Vogler et al., 2012;Tusso et al., 2016). Thus, this gene marker may be too conservative for detecting haplotype and nucleotide diversity, which likely led to the results of the genetic diversity analysis of COTS based on the mtDNA gene marker not being consistent with those of the microsatellite loci analysis. However, the mtDNA marker could be used to analyze the demographic history and lineage divergence of COTS populations (Vogler et al., 2008(Vogler et al., , 2013. None of the COTS populations in the SCS showed lineage divergence from that in Pacific Ocean, which may have also diverged in the Pliocene-Early Pleistocene (1.95-3.65 Mya; Vogler et al., 2008). The speciation process of Pacific, Red Sea, and Southern and Northern Indian Ocean COTS groups was likely driven by largely biogeographic barriers between the main oceans due to sea-level changes (Pillans et al., 1998). Moreover, Tajima's D, Fu' Fs, and mismatch distribution indicated that COTS experienced historical population expansion in four coral habitats in the SCS (Table 1 and Figure 2C), which was consistent with the demographic history of COTS populations in the Western Pacific (Vogler et al., 2013). These results may be closely associated with coral population expansion that was driven by sea-level rise after the Last Glacial Maximum (Voris, 2000;Tian et al., 2010;Huang et al., 2018). In the Pleistocene, glacial advances led to a sea-level fall to 120-140 m below present levels in Southeast Asia and Australasia, which led to the bulk of the Sunda and Sahul shelves being largely exposed and formed massive lowland connections between present day islands in this region and adjacent continents (Voris, 2000). Meanwhile, coral habitats in the northern SCS also disappeared during this period (e.g., Weizhou Island, Hainan Island, and Daya Bay) due to land advances (Voris, 2000;Huang et al., 2018), and many marine organisms, including corals, retreated into periodic refugia and subsequently expanded after the end of the glacial age. For example, Huang et al. (2018) found that populations of Porites lutea expanded after the Last Glacial Maximum in the SCS, and re-established populations on Weizhou Island and Hainan Island from the northern SCS as the sea levels rose.
COTS live in coral habitats and feed on scleractinian corals (Glynn, 1988;Rotjan and Lewis, 2008), and consume a maximum of 66-478 cm 2 of live coral per day, depending on season and their body size (Glynn, 1973;Keesing and Lucas, 1992). However, densities of 10-15 starfish per hectare (0.001-0.0015 ind./m 2 ) could be sustained in areas with more than 20% coral cover (Keesing and Lucas, 1992). Therefore, re-established coral may provide more food sources for COTS in a larger biogeographical range, and it is reasonable to speculate that the COTS population historical expansion followed coral habitat range expansion. Moreover, COTS generally prefer to prey on Acropora spp. and Pocillopora spp. (Kayal et al., 2012), which were widely distributed in the mature tropical coral reef ecosystems in the SCS, and were rare inshore (Chen et al., 2009;Zhao et al., 2013Zhao et al., , 2016Yu et al., 2019). Birkeland and Lucas (1990) found that variation in the effects of COTS outbreaks in distinct coral habitats of the Pacific may be explained by the relative dominance of Acropora in coral communities. Furthermore, Acropora spp. are consistent among the corals that are worst affected during outbreaks and tend to dominate in the Western and Southern 3 | Analysis of molecular variance (AMOVA) for crown-of-thorns starfish in the South China Sea and Pacific Ocean from each 17 coral habitats with 4 geographical regional populations (South China Sea, Western Pacific, Pacific equatorial currents affected zone and Pacific insular atoll).  Pacific (Pratchett et al., 2009;Pratchett, 2010). Thus, feeding preferences of COTS may lead to COTS population expansion relying on coral population range expansion. It was found that, as yet, there were no spines or individual COTS in subtropical coral communities in the northern SCS, which are dominated by massive corals, such as Dipsastrea palauensis, Porites lutea, and Dipsastrea halicora (Chen et al., 2009;Liao et al., 2019;Yu et al., 2019). Furthermore, the large exposure of Sunda and Sahul shelves and the decrease in width of the Balabac Straits (∼12 km) reduced marine organism connectivity between the SCS and surrounding oceans in the glacial period of the Pleistocene (Voris, 2000). However, the Luzon Strait was not closed due to land expansion and was always connected between the SCS and Western Pacific (Voris, 2000), which was also the main seawater exchange channel in the SCS (Yu et al., 2010). Thus, the opening of the Luzon Strait may be the main reason for there being no genetic divergence in the conserved mitochondrial COI fragment of COTS populations between the SCS and the Pacific.
Accordingly, COTS may have undergone historical population expansion after the Last Glacial Maximum in the Pleistocene, which might have been followed by the expansion of coral populations in the SCS. Furthermore, the biogeographical range of COTS population expansion may have been limited by feeding preferences.

Genetic Differentiation of Crown-of-Thorns Starfish Populations May Be Shaped by Distinct Nutrient Concentration in the South China Sea
It is difficult to determine the relationship between environmental factors and the genetic structure of COTS, as COTS cannot be accurately aged (Deaker et al., 2020). However, this problem can be partially solved by measuring the diameter of COTS (Moran, 1988;Pratchett et al., 2014) and selecting an appropriate sampling time (such as a population outbreak period; Pratchett et al., 2009;Saponari et al., 2018). The nutrient concentration analysis results showed that the Chl a concentration in DI was significantly higher than that in the other coral habitats in the SCS (Figure 6A). Additionally, the POC concentration in DI and BR was significantly higher than that in FTS and PK (Figure 6B), which was consistent with the genetic clustering results inferred by STRUCTURE ( Figure 4B). Furthermore, the Mantel test indicated that there was a significant correlation between the Fst values and POC concentration in the SCS (Figure 5C), which was consistent with the DistLM analysis results. Thus, the significant genetic differentiation of COTS in distinct coral habitats with slight geographic isolation may be shaped by nutrient concentration, particularly the POC concentration. One adult female COTS can produce approximately 10 million fertilized eggs per spawning (Conand, 1984;Babcock et al., 2016;Caballes et al., 2021;Pratchett et al., 2021), and COTS larvae predominantly feed on mid-sized phytoplankton (e.g., dinoflagellates and pennate diatoms > 5 µm; Ayukai, 1994). However, the larval starvation hypothesis suggests that severe food limitation would lead to high mortality in COTS larvae (Lucas, 1982;Fabricius et al., 2010). Thus, distinct nutrient concentrations may cause stress on COTS larvae among coral habitats, which may lead to genetic differentiation. Pratchett et al. (2017a) combined many studies and suggested that Chl a concentrations between 0.4 and 1.0 µg·L −1 may be optimal for COTS larval survival and settlement. Nevertheless, the total food abundance and sources for COTS larvae still cannot be accurately indicated by Chl a concentrations (Ayukai, 1994;Tada, 2003;Mellin et al., 2017;Pratchett et al., 2017a), COTS larvae not only selectively ingest some phytoplankton species (e.g., Chaetoceros sp. and Dunaliella sp.; Ayukai, 1994;Mellin et al., 2017), but can also utilize other nutrients in oligotrophic coral reefs (Olson and Olson, 1989;Hoegh-Guldberg, 1994;Nakajima et al., 2016). These were the main reasons for the lack of a significant correlation between the Chl a concentration and genetic differentiation indices (Fst) of COTS in the SCS.
It should be noted that COTS larvae can utilize nutrients directly associated with POC and dissolved organic carbon. Some studies found that COTS larvae can assimilate and utilize dissolved organic matter and bacteria in the sea water column (Olson and Olson, 1989;Hoegh-Guldberg, 1994), and the large amounts of energy required for developing larvae could be provided by dissolved free amino acids (Hoegh-Guldberg, 1994). However, the concentration of dissolved free amino acids in natural seawater is extremely low, and cannot make a significant contribution to the nutrition of COTS larvae (Ayukai, 1993). Moreover, coral and other reef benthos can release large amounts of mucus aggregates and organism-derived organic matter to reef water, which are closely associated with POC and provide alternative food resource for COTS larvae. Many previous studies also showed that particulate organic matter in reef water can be utilized by COTS larvae (Okaji et al., 1997). An experiment with labeled stable isotope tracers ( 13 C and 15 N) also found that COTS larvae can utilize coral mucus particles and associated microorganisms, which may be an additional food resource (Nakajima et al., 2016). Therefore, the impact of Chl a concentration on the genetic structure of COTS larvae in the oligotrophic waters around Paracel and Spratly Islands (Nansha Islands) was weak. But, distinct POC concentration may directly lead to significant genetic differentiations of COTS populations among coral habitats in the SCS, over small geographic distances (e.g., DI-BR genotype and FTS-PK genotype; Figure 4). In addition, because the spatial difference of nutrient availability has a substantial impact on the survival and growth of COTS larvae (Olson, 1987;Okaji et al., 1997;Wolfe et al., 2015), it was speculated that COTS larvae of the FTS-PK genetic cluster may be able to tolerate low nutrient concentrations, whereas DI-BR COTS genotypes are more competitive in coral habitats with higher POC concentration. The low nutrition resistance of COTS larvae may be associated with maternal nutrition . Previous studies have shown that the coral cover of Acropora was decreased with the increasing of the latitude in the SCS, and the abundance of Acropora in FTS and PS was higher than that in DI and BR (Zhao et al., 2013(Zhao et al., , 2016Li et al., 2019;Liao et al., 2021). Acropora was as preferred coral prey for COTS, which may provide high quality and quantity of food to female COTS in the FTS and PS, which might promote COTS larvae to tolerate low nutrient availability and to progress to more advanced stages faster . However, this advantage may diminish under the high POC condition.
There was a weak genetic differentiation in COTS populations between FTS and PK, which implied that COTS larval populations have effective connectivity between two coral habitats with a distance of ∼800 km (Figure 1). In addition, the STRUCTURE results showed that the frequency of the FTS-PK genotype (blue) was also higher than that of the DI-BR (gray) genotype in the COTS populations in the Western Pacific. Thus, the FTS-PK genotype may have a longer PLS, which might be able to achieve long-distance dispersal via ocean currents. It has been found that several COTS larval genotypes have a relatively long PLS period in the Pacific Ocean, and have spread among Guam, Kingman Reef, and Swains through equatorial current systems (Tusso et al., 2016). However, the frequency of the DI-BR genotype was rare except on the Paracel Islands, which suggested that the DI-BR genotype may belong to the local COTS population.
Accordingly, distinct nutrient concentration, especially different POC concentration may lead to genetic differentiation of COTS populations among coral habitats within short geographic distance in the SCS. Also, FTS-PK genotypes may have relative longer PLS and potentially stronger dispersal ability than the DI-BR genotype in the Paracel Islands.

Ocean Currents Provide a Potential Dispersal Mechanism for Crown-of-Thorns Starfish Larvae in the South China Sea
The COTS populations have effective connectivity and gene flow between PK and FTS (Figure 4B), which may benefit from the transport of COTS larvae by ocean currents. The COTS are similar to many reef benthic organisms (e.g., coral and anemones; Kool et al., 2011;Treml et al., 2012), which are almost incapable of migrating between isolated atolls as adults (move < 35 m/day; Chesher, 1969;Keesing and Lucas, 1992;Pratchett et al., 2017c). Thus, the dispersal of COTS was more likely to be completed in the PLS. The shortest geographic distance is about 800 km between PK and FTS, which is a large geographic barrier for COTS larval dispersal. The PLS of COTS lasts between 9 and 42 days , and therefore, larvae were expected to disperse only 10-100 km between reefs without assistance from oceanic currents (Dight et al., 1990). The results of IBD analysis also suggested that large geographic isolation would lead to a genetic break among COTS populations ( Figure 5A; Dight et al., 1990;Vogler et al., 2012). For instance, there was significant genetic differentiation between COTS in Southeast Africa and that in East Asia (Yasuda et al., 2009). However, the gene flow inference suggested that COTS populations have effective connectivity between reefs separated by <1000 km (Yasuda et al., 2009). Unexpectedly, it has been found that COTS larvae can spread to remote coral reefs with the assistance of strong ocean currents in the Pacific Ocean (Vogler et al., 2013;Tusso et al., 2016). In addition, the spawning period of COTS is from May to August in the Northern Hemisphere (Pratchett et al., , 2017a, and thus, the ocean currents from May to August were critical to COTS individual dispersal in the SCS. The results of ocean current analysis showed that there were stable currents from the southern part of the Paracel Islands to the Spratly Islands from June to August (current velocity: 0.25-0.4 m/s; average current velocity: 0.33 m/s; Figure 7), which may assist COTS larvae in the Paracel Islands waters to travel to the Spratly Islands within 23 days. Therefore, based on these results, we suggest that some COTS larvae in FTS may have originated from the southern part of the Paracel Islands. Since the ocean current was from south to north in May from the SCS, the migration of COTS larvae from PK to FTS may have started in June.
Interestingly, the STRUCTURE results indicated that blue and green genotypes were shared between the FTS-PK and Western Pacific, which implied that COTS populations exhibit recombination and gene exchange between the two regions ( Figure 4B). The SCS is a semi-enclosed marginal sea, and the main seawater exchange is through the Luzon Strait to the Pacific Ocean (∼2400 m; Yu et al., 2010). The Taiwan Strait (∼70 m) and the Mindoro Strait (∼400 m) are outflow zones for the SCS (Fang et al., 2003;Su, 2004;Yu et al., 2010). Since the Karimata Strait is shallow (∼29 m), there is negligible exchange of water between the SCS and the Indian Ocean (Yu et al., 2010). Some studies have suggested that reef organisms in the southern SCS have genetic connectivity with the "Coral Triangle, " and Symbiodiniaceae and coral larvae could enter the SCS through the Sulu Sea (Kool et al., 2011;Chen et al., 2020). However, many marine organisms that disperse mainly during the PLS still have large biogeographical barriers between the SCS and the Indo-Pacific region. Notably, recent gene flow and ocean current data revealed that western boundary currents and coastal currents play an important role in COTS larval dispersal in the Pacific Ocean (Dight et al., 1990;Moran et al., 1992;Yasuda et al., 2009). For instance, COTS larvae were transported from Mindoro in central Philippines to northern Luzon Island by the Northwest Luzon Coastal Current (Hu et al., 2000;Yasuda et al., 2009), and the strong Kuroshio Current is likely to continue to transport them to Ryukyu Island in southern Japan (Yasuda et al., 2009;Yasuda, 2018). Some mtDNA and microsatellite data have shown that there was high connectivity of COTS populations between the Philippines and the Western Pacific Ocean (Vogler et al., 2013;Tusso et al., 2016). Thus, the boundary and coastal current systems may help the COTS larval populations in the Western Pacific enter the SCS. Many COTS larvae may have migrated and spread among many coral habitats in a stepping-stone fashion within the SCS (Tusso et al., 2016).

CONCLUSION
Based on mtDNA sequences, and phylogenetic and demographic history analysis of COTS populations, there was no significant phylogenetic divergence in COTS populations between the SCS and Pacific Ocean. However, there was a clear signature of population expansion in the SCS. These results suggested that COTS may have undergone historical population expansion after the Last Glacial Maximum in the Pleistocene, which might have been followed by the expansion of coral populations in the SCS. Furthermore, the opening of the Luzon Strait may be one of the reasons for the lack of significant genetic divergence in COTS populations between the SCS and the Pacific. In addition, according to the microsatellite loci analysis, the genetic diversity of COTS populations in the SCS was consistent with that in the Western Pacific and Pacific equatorial current affected zone COTS populations, which was higher than that in the Pacific insular atolls. However, the genetic structure of COTS populations in the SCS may have been shaped by distinct nutrient concentrations, particularly different POC concentrations over small geographic distances. FTS-PK genotypes may have stronger dispersal ability than the local DI-BR genotype in the SCS. Moreover, ocean currents provide a potential dispersal mechanism for COTS larvae in the SCS. Stable ocean currents may assist COTS larvae to travel from the Paracel Islands to the Spratly Islands from June to August, and COTS larvae in the Western Pacific may undergo gene exchange with those in the SCS by boundary and coastal current systems. Accordingly, this study demonstrates that environmental and oceanographic factors play important roles in shaping the genetic characteristics and larval dispersal of COTS populations in the northern Indo-Pacific convergence region.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/ Supplementary Material.

AUTHOR CONTRIBUTIONS
KY and BC deigned the research. QY, ZL, ZQ, XY, QW, and BH contributed the materials. BC, QW, and BH performed the research. BC analyzed the data and drawn all pictures. BC and KY wrote the manuscript. All the authors reviewed the manuscript.

FUNDING
This work was supported by the National Natural Science Foundation of China (Nos. 42090041 and 42030502), the Guangxi Scientific Projects (Nos. AD17129063 and AA17204074), and the Innovation Project of Guangxi Graduate Education (No. YCBZ2018006).