The potential role of temperate Japanese regions as refugia for the coral Acropora hyacinthus in the face of climate change

As corals in tropical regions are threatened by increasing water temperatures, poleward range expansion of reef-building corals has been observed, and temperate regions are expected to serve as refugia in the face of climate change. To elucidate the important indicators of the sustainability of coral populations, we examined the genetic diversity and connectivity of the common reef-building coral Acropora hyacinthus along the Kuroshio Current, including recently expanded (<50 years) populations. Among the three cryptic lineages found, only one was distributed in temperate regions, which could indicate the presence of Kuroshio-associated larval dispersal barriers between temperate and subtropical regions, as shown by oceanographic simulations as well as differences in environmental factors. The level of genetic diversity gradually decreased towards the edge of the species distribution. This study provides an example of the reduced genetic diversity in recently expanded marginal populations, thus indicating the possible vulnerability of these populations to environmental changes. This finding underpins the importance of assessing the genetic diversity of newly colonized populations associated with climate change for conservation purposes. In addition, this study highlights the importance of pre-existing temperate regions as coral refugia, which has been rather underappreciated in local coastal management.

are classified as Near Threatened or Vulnerable on the IUCN Red List, temperate regions may serve as refugia for hermatypic coral species in the near future 3 if the coral species can overcome the seasonal fluctuations in temperature and reduced light intensities 5 . Among the multiple definitions of refugia, we use the term as defined by Keppel and Wardell-Johnson 6 . In this paper, authors noted that "Refugia are habitats that components of biodiversity retreat to, persist in and can potentially expand from under changing environmental conditions and facilitates the persistence of species under anthropogenic climate change. " Population connectivity and relatively high genetic diversity are crucial for a population to expand and persist under climate change 7 ; thus, these characteristics are important for the populations in potential refugia.
Winter water temperature has risen by 1.45 °C over the last 100 years in the temperate Kuroshio Current region in Japan (http://www.data.jma.go.jp/gmd/kaiyou/data/shindan/a_1/japan_warm/japan_warm.html), causing coral expansion and producing potential coral refugia at the edges of their distribution ranges 3 . While such temperate region can be potential coral refugia, temperate peripheral populations may suffer from reduced genetic diversity due to small population sizes and isolation 8 . However, no studies have assessed the genetic diversity or structures of such coral populations, including the recently colonized temperate populations that expanded within the last 50 years as a result of climate change along the Kuroshio Current. In addition, there is a lack of knowledge of the oceanographic larval dispersal processes that shape the genetic variations and population structures of corals in the high-latitude regions of Japan.
This study aimed to assess the genetic structures of recently expanded populations and examine the potential role of temperate regions as coral refugia in Japan. We thus examined the genetic diversity and meta-population structures of the subtropical to northernmost temperate populations of the reef-building coral Acropora hyacinthus, and the results were further compared with contemporary larval dispersal simulations in an oceanographic model.
Our target species, Acropora hyacinthus (Dana, 1846), is widely distributed in the tropical, subtropical, and some part of temperate regions of the Indo-Pacific 9 . Northward range expansion of A. hyacinthus has been observed in the temperate regions of Japan for the last 50 years 3 . Previous studies that included some temperate populations found cryptic lineages (species) of A. hyacinthus 10,11 . The A. hyacinthus species complex is a hermaphroditic, broadcast-spawning species whose fertilization takes place in the water column. More than 70% of the larvae settle within one week 12 but can survive in the water column for up to 5-6 weeks 13 . A. hyacinthus is categorized as Near Threatened on the IUCN red list. A. hyacinthus is one of the key reef-building species on tropical reefs 14 and can represent as much as 36% of the coral coverage in some high-latitude areas 15 . These ecological features of this species are suitable for comparing the genetic diversity and connectivity of populations among temperate and subtropical regions. We defined the region south of the Watase Line 16 as subtropical and the region north of the Watase Line as temperate (Fig. 1). The region was split along this line because the biogeography of this region is often split between Yakushima and the Amami Islands, and the Kuroshio Current passes through this region.

Results
the Acropora hyacinthus lineage used in this study. Only 10366m5 has possible null alleles, and we thus excluded this locus from the genetic analysis of the population. A previous study indicated that there are cryptic species (lineages) in Acropora hyacinthus that are morphologically indistinguishable 10 . Thus, we identified the cryptic lineages in our dataset using STRUCTURE analysis 17 to determine the lineage associated with the recently expanded populations in temperate regions. After removing 44 possible clones, our STRUCTURE analyses based on 645 multilocus genotypes (MLGs) indicated that the most likely number of clusters was K = 2, followed by 3 (Suppl. 1). At K = 2, temperate populations and some of the subtropical individuals were assigned to the same cluster (green lineage in Fig. 2), and the other lineage was restricted to south of Miyazaki. At K = 3, the subtropical specific cluster at K = 2 (red in Fig. 2) was further split into two different clusters (red and blue in Fig. 2). The principal coordinate analysis (PCoA) results also supported different lineages (Fig. 2). Moreover, populations with two or three different genetic lineages belonged to different clades in the neighbour-joining (NJ) phylogenetic trees (Suppl. 2). These results, as well as those of a previous study 11 , indicated that these populations are possibly different species 10 . Because we are interested in the genetic connectivity and diversity of the lineages that have recently expanded in temperate regions, we focused on one lineage (green lineage in Fig. 2) found in both temperate and subtropical regions. Individuals assigned to the green lineage with probabilities ≥70% in both the K = 2 and K = 3 results were used for the following analysis (Table 1). A relatively low value ≥70% was chosen because that lineage was rather rare in the subtropical region and many individuals found in the subtropical area are partially admixed with other lineages. Populations with ≥9 MLGs were used for subsequent analysis because smaller sample sizes are not suitable for population genetic analysis 18 . The final dataset included 240 individuals from 13 populations (shown in bold in Table 1). Initially, significant linkage disequilibrium was found in two pairs of loci (11401m4 and Amil2-23; 11401m4 and Amil2-2). However, linkages disappeared in the final dataset, which included only one lineage and no clones.
Genetic diversity. Possible clones with a probability of identity of P < 0.001 were more frequently found in the temperate region than in the subtropical region (9 out of the 11 populations in the temperate region but none in the subtropical region), which indicated that asexual reproduction (such as fragmentation) likely occurred more frequently in the temperate region than in the subtropical region (Table 1). High clonal rates (23.1% for GAHk and 52.6% for GAHo) were found in two Goto Island populations, both of which are recently expanded populations.
After removing possible clones, the north peripheral populations had reduced genetic diversity, which was observed in particular in the recently expanded populations (Fig. 3, Suppl. 3). Genetic diversity metrics (allelic richness and heterozygosity) were more significantly reduced in the recently expanded populations than in the rest of the pre-existing temperate populations (Wilcoxon test, P = 0.010 and P = 0.038 for allelic richness and expected heterozygosity, respectively) (Fig. 3). Similarly, the genetic diversity metrics were significantly higher in subtropical populations than temperate ones (Wilcoxon test, P = 0.007 and P = 0.028 for allelic richness and expected heterozygosity, respectively). Notably, all temperate populations, except for one north Yakushima population (KAHyn), had more than one fixed (monomorphic) locus. All the recently expanded populations had 2 to 4 fixed loci out of 8 loci (Suppl. 3). In total, 24 individuals had private alleles. Eight out of the 32 individuals (25.0%) examined in the subtropical populations had private alleles, and 11 out of the 149 (7.4%) individuals in the pre-existing temperate populations had private alleles. Two out of 59 individuals (3.4%) were found to have private alleles in the recently expanded populations (One from Goto Island, GAHo, and another from Amakusa, Isolation by distance (IBD) pattern and genetic connectivity. Significant IBD patterns within the green lineage were detected only when regressed onto oceanographic distance (P < 0.05) (Fig. 4, Suppl. 4). A weak IBD tendency was detected using straight-line distance (P = 0.085) only when the recently expanded populations were excluded, which possibly occurred because the recently expanded populations are not in genetic equilibrium. The analysis of molecular variance (AMOVA) results suggested overall significant genetic structuring within the A. hyacinthus lineage that was examined in this study (global F ST = 0.108, P < 0.001, Table 2). Significant genetic differentiation between the temperate and subtropical populations was not found, either including (all temperate vs subtropical, F CT = 0.042, P = 0.064) or excluding (pre-existing temperate vs subtropical, F CT = 0.024, P = 0.163) the recently expanded populations. The AMOVA results revealed non-significant   (Table 2). This result suggests that populations in temperate regions were generally genetically close to each other. Pairwise tests showed almost no genetic differentiation among the recently expanded populations, whereas there were some significantly differentiated pairs among the pre-existing populations ( Table 3). The Tsukijima population in Miyazaki (TUK) was one of the most isolated populations in the temperate region, most likely because the artificial structure surrounding this population (N. Yasuda pers. obs.) hindered larval exchange. STRUCTURE analysis was again performed within the green lineage with and without prior information (subtropical vs temperate). The results indicated that K = 1 best explained the data (Suppl. 5a,b), which suggests admixture and no clear genetic divergence between the temperate and subtropical regions, and these findings are consistent with the AMOVA results. However, when using the LOCPRIOR model, different proportions of assigned lineages were found between temperate and subtropical regions at K = 2 (Fig. 5) with a small r value (0.869), which implies that very small genetic differences exist between the two regions.
The estimation of recent gene flow based on the assignment method using GeneClass2 showed that 95 out of 240 individuals were successfully assigned to one of the populations. Most of the successfully assigned individuals (15 out of 16) in the recently expanded populations suggested self-seeding. Only one out of 76 (1.3%) temperate individuals were possible recent immigrants from subtropical populations (Suppl. 6). particle tracking simulation. Larval dispersal generally occurred northward along the Kuroshio Current (Fig. 6b,c). Larval dispersal from temperate to subtropical areas seemed uncommon within a single generation because larval dispersal is physically limited to the path of the Kuroshio Current (here called the Kuroshio-associated barrier), although dispersal could be achieved through multiple generations or occasional  long dispersals (Fig. 6c). The Kuroshio-associated barrier corresponds to the biogeographic boundary called the Watase Line (Fig. 1). Unidirectional dispersal was observed from the western Kyushu region (E in Fig. 6a) along the Tsushima Current to the temperate region along the southern Japanese coast (B in Fig. 6a) along the main stream of the Kuroshio Current (Fig. 6c).

Discussion
Cryptic lineages were found in A. hyacinthus, while only one lineage is dominant in the temperate region and is apparently expanding to the north. The Kuroshio-associated barrier found in the numerical simulation and/or the lack of adaptability to the temperate environment could hinder northward migration of the blue/red lineages (Fig. 2). The Kuroshio-associated barrier also matched the Watase Line, which is the biogeographic barrier that separates terrestrial fauna and was caused by land bridges and sea barriers in past climates 16 . Larval dispersal simulations using particle tracking averaged over 24 years indicated discontinuous patterns of connectivity between Yakushima and the Amami Islands (Fig. 6). This discontinuity is possibly due to the Kuroshio path crossing between Yakushima and the Amami Islands, which prevents direct dispersal of single generations from the subtropics to areas further north (referred to as the Kuroshio-associated dispersal barrier). This study demonstrates that the Watase Line can be found in marine species due to the existence of the Kuroshio-associated dispersal barrier. The limited larval transport from subtropical to temperate regions within a single generation (Fig. 6) suggested that the probability of future poleward migration of the other lineages (red/blue in Fig. 2) of A. hyacinthus would depend on whether the larvae of these lineages could pass through this barrier via long or multi-generation   dispersal. Further regular monitoring of the distributions of all cryptic lineages using genetic markers is necessary to assess the poleward migration of the species. Reduced genetic diversity (heterozygosity, allelic richness, and clonal rate) was observed in the recently expanded Acropora hyacinthus populations distributed near the edge of the species' range. One historical interpretation is that reduced genetic diversity near the edge of a species' range is the result of greater isolation in peripheral populations than in core populations where species are more abundant than along the periphery 8 . However, this study did not detect a gradual increase in genetic isolation towards peripheral populations (Table 3). Alternatively, IBD patterns produced by a simple linear stepping-stone model explains the reduced genetic diversity near the edge. If a population structure is assumed to exhibit linearly uniform distribution over a large distance with limited migration in a simple linear stepping-stone model, then IBD genetic structure can theoretically be observed 19 . Under such a situation, Tajima 19 demonstrated that the level of DNA polymorphism is highest in the central populations and gradually decreases towards peripheral populations. Our data relatively fit this linear stepping-stone model because significant IBD patterns were observed when oceanographic distance was used. Although more data are needed, a linear stepping-stone model can at least partially explain the observed patterns. Lastly, it is also possible that genetic diversity will decrease towards the edge of a species' range in response to different selection pressures 20 , although there is an opposing view that a stressful peripheral environment can increase the number of genotypes in some cases 21 . Environmental factors, such as low sea water temperature 22 , high concentration of chlorophyll a, particulate organic carbon (POC) 23 , reduced aragonite saturation 24 , light intensity 25 , and high turbidity 26 around northern populations (Suppl. 7), could affect the suitability of reef-building corals habitats 25 and subsequently erode genetic diversity. On the other hand, the larger genetic diversity in the subtropical populations can be partly explained by possible hybridization between closely related lineages (as found in Fig. 2). Some previous studies have also provided evidence for reduced genetic diversity in peripheral coral populations. For example, both brooder and broadcaster coral species in the southernmost reef in the Pacific showed lowered genetic diversity 27 . In the Atlantic, another coral species, Montastraea cavernosa, showed lower genetic diversity in peripheral populations 28 . These studies suggested that a small founding population, small effective population size of genetic drift, or inbreeding depression caused reduced genetic diversity in peripheral coral populations. However, expanding coral populations do not always follow this trend of reduced genetic diversity, as was demonstrated in the Mediterranean coral Oculina patagonica 29 .
The clonal rates were higher in temperate regions than in subtropical regions. Although the mechanism is unclear, different habitat characteristics (i.e., reef structure in subtropical regions and non-reef structure in temperate regions) potentially caused different rates of sexual reproduction success in A. hyacinthus. In Acropora palmata, different rates of asexual reproduction are known to be related to different habitat characteristics 30 .
While the historical gene flow analysis (e.g., F-statistics; Tables 2 and 3, Suppl. 5) showed overall genetic similarities among temperate populations, contemporary gene flow analysis (assignment tests) indicated that all studied populations, including those in the peripheral temperate regions, are basically self-seeding (Suppl. 6). This result concords with the fact that peripheral temperate coral populations are mostly self-sustained by local spawning 31 . Recently expanded populations were therefore possibly seeded from other temperate populations by occasional long dispersal and then sustained themselves by self-seeding.
Although the recently expanded Goto populations (GAHo) had low genetic diversity and a small number of samples was analysed (MGL = 9), one individual had a private allele. One possible explanation is that harsh environments around GAH (Suppl. 7) caused some strong directional selective pressures on the coral populations, which reduced the genetic diversity and changed the allele frequencies, as was discussed above (selective pressure increased the frequency of rare alleles, and a private allele was found). Alternatively, the Goto population was colonized by other genetically diverged populations, such as those from Taiwan through the Taiwan-Tsushima warm current system 32 , that were not analysed in this study. Previous study has claimed that peripheral populations with low genetic diversity are some of the most active regions of speciation because of such genetically distinct features 33 . Further genetic analysis of Goto will be interesting to examine potential speciation in such a peripheral region.
Significant IBD patterns were detected only when oceanographic distance was used (Fig. 4), which indicates that pairwise genetic distances (F ST ) significantly increased with geographic distance. This result indicates that overall, the effect of genetic drift is stronger than the effect of gene flow between distant populations.
We expected that the larval dispersal model rather than oceanographic distance would have significant IBD patterns because genetic structure is considered to be most associated with larval dispersal. However, the IBD patterns that were derived using the number of larvae between sites based on numerical simulations did not improve in comparison to those using Euclidean distances and oceanographic distances (Fig. 4, Suppl. 4). This result might be because the genetic markers used in this study captured historical events rather than contemporary gene flow, especially in the F ST calculations, which assume equilibrium. Similarly, the genetic distance of sea stars across the Indo-Pacific region was best explained by Euclidean distance and not modelled larval distance in a previous study 34 in which the authors considered that the estimated F ST reflected historical gene flow rather than contemporary gene flow. Another explanation for the weak correlation with IBD patterns is the lack of biological characteristics (i.e., larval behaviour) in the larval dispersal models used in this study. A more sophisticated model that includes the biological aspects of A. hyacinthus would improve the correlation with genetic distance, as in a previous study 35 .
Although corals in subtropical regions are threatened due to climate change, they are expanding and propagating into temperate regions 3 . Many climate change predictions of the species distributions of corals show a general loss of reef habitats in many areas except for the relatively robust eastern marginal Pacific region, including temperate Japan 36 . Given the moderate levels of neutral genetic diversity and the persistence and increase in coral cover before and after the rise in sea water temperatures in pre-existing temperate coral populations, we propose that some of the pre-existing temperate regions in Japan would be good refugia for at least one of the A. hyacinthus lineages. Pre-existing temperate regions are indeed the potential habitats that "part of biodiversity retreat to, persist in and can potentially expand" under climate change conditions. In temperate regions, seaweed beds have been replaced by temperate corals, including A. hyacinthus, and the coral cover is currently expanding in these regions 37 . Although a direct association between neutral genetic diversity and resilience of coral populations against environmental change is partly under debate 8 , some coral genotypes are known to be associated with stress tolerances such as disease and higher water temperature 38 . Thus, moderate genetic diversity in pre-existing temperate populations may indicate stability and higher tolerance against unexpected extreme environmental changes such as local anthropogenic stress. Nonetheless, much less attention has been paid to coral assemblages and their conservation in temperate regions than in subtropical regions. Unlike in subtropical regions, where tourism, recreation, and fisheries are directly associated with the abundance of corals, temperate corals are sometimes regarded as a nuisance by local fishermen in temperate regions because subtropical fishes are of low economic value in these areas. In addition, expensive fishing nets become entangled on and broken by hard corals, causing pecuniary damage to fishermen (N. Yasuda personal communication with local fishermen in Miyazaki). Because reef-building corals help to increase overall biodiversity 39 by providing food and habitats with complex structures in coastal ecosystems, there is a need to review the value of temperate corals for conservation.
The importance of conserving the pre-existing temperate corals discussed in this study does not dismiss or discourage the conservation of subtropical corals. Rather, two out of the three A. hyacinthus lineages were found in only the subtropical region, and the highest genetic diversity was found in the subtropical region, indicating that these areas are nonetheless important.
This study indicated that recently expanded marginal populations could have low genetic diversity and are possibly more vulnerable to environmental changes than pre-existing populations. This finding underpins the importance of assessing the genetic diversity of populations that have recently expanded due to climate change for conservation purposes, in particular when anticipating the roles of these areas as refugia under climate change.

Methods
Coral sampling and genotyping. In total, 689A. hyacinthus samples were collected using SCUBA from 22 populations along the Kuroshio and Tsushima Currents (Table 1, Fig. 1). This study included four temperate populations that have recently expanded possibly due to climate change: two populations from Goto Island and the Amakusa and Shikine populations 3 . No records of A. hyacinthus were found in Goto and Amakusa in a 1965-1977 survey. There was no existing record of A. hyacinthus in Shikine in a 1986 survey 3 , but we discovered large colonies in 2016 in this study. Each sample was photographed except for the AMK, KOC, MYK, SSK, TUK, TOI, and OAH specimens (Table 1). Wherever possible about 5 cm of skeleton specimens were collected and deposited at the University of Miyazaki (MUFS-Acr101-310, 371-428, 500-531, 569-717, 768-801). Species identification was further confirmed by Dr. Hironobu Fukami based on morphological information. Genomic DNA was extracted using a modification of a previously published method 40 . Genotypes of the corals were determined by 8 microsatellite loci 41,42 . Eight loci were amplified in 2 multiplex PCR sets using different dyes and non-overlapping loci (Plex1 includes 8346m3, 11401m4, 8499m4, and 10366m5; Plex2 included Amil2-2h, Amil2-23h, Amil2-10, and Amil2-12). For Amil2-2h and Amil2-23h, we used modified primers to increase PCR efficiency for A. hyacinthus (Suppl. 8). We used an 8 μl reaction containing 1 μl of template DNA, 0.03 μl of each primer (50 μM), 4 μl of Type It (Microsatellite Kit, Qiagen), and 2.94 μl of water. The thermal cycler profile that was used had an initial denaturation step at 95 °C for 3 min followed by 30 cycles of 95 °C/30 s, 50 °C/30 s, 72 °C/30 s, with a final extension step at 72 °C for 5 min.
The multiplexed PCR product was diluted by 100, and 1 μl of the diluted PCR product was added to 10 μl of highly deionized formamide containing 0.15 μl of the GeneScan-600 (LIZ) size standard (Applied Biosystems, Foster City, CA, USA) for genotyping on an ABI 3730xl automated sequencer (GeneMapper ® ver. 5, Applied Biosystems), which was used to determine the fragment sizes (alleles) of all samples. All genotyping was manually checked. Whenever we found ambiguous patterns or no peak, we re-extracted the genomic DNA and re-amplified it or removed the sample from the analysis. We tested for the presence of null alleles using MICRO-CHECKER ver. 2.2.3 43 , and the locus positive for null alleles was excluded from the genetic analysis of the population.
The clonal structure of each coral population was examined by identifying the same MLGs based on 8 microsatellite loci; these were identified as possible clones (asexually reproduced by fragmentation), and all except one were excluded from subsequent analyses. The inbreeding coefficient of a population (F IS ), was evaluated with GenAlEx ver. 6.502 44 . Genotypic linkage disequilibrium across 8 loci was examined using GENEPOP on the web (ver. 4.0.10) 45 by estimating exact P values using the Markov chain method with 10,000 dememorization steps, 1,000 batches, and 10,000 iterations per batch.
Different measures of genetic diversity were estimated for each population, including the number of alleles per locus (Na), the observed (H O ) and expected (H E ) heterozygosity using GenAlEx ver. 6.502 44 , and allelic richness (A) using FSTAT ver. 2.932 46 .
Identifying cryptic lineages distributed in both temperate and subtropical regions. We first identified the cryptic lineages of A. hyacinthus found in a previous study 11 . STRUCTURE analysis was conducted using the admixture model without prior information. The Markov chain Mote Carlo (MCMC) chains used a burn-in of 500,000 chains followed by 500,000 MCMC replications. Ten independent chains were run for each K = 1 to 10 to confirm convergence. The coefficient of ancestry was calculated and visualized for each individual across 10 runs for the most likely values of K (Delta K 47 ), in CLUMPAK 48 Among the lineages identified in the STRUCTURE analysis above, we used only the lineage distributed in both temperate and subtropical regions in subsequent analyses (green lineage in Fig. 2). PCoA was conducted to confirm the different lineages found in the STRUCTURE analysis using GenAlEx ver. 6.502 44 . A NJ phylogenetic tree was reconstructed using Nei's genetic distance 49 as implemented in POPULATIONS ver. 1.2.32 50 with 1000 bootstrap values and visualized this tree using TreeViewX ver. 0.5.0 51 . STRUCTURE was again run to examine the intra-green-lineage genetic structure. First, we analysed STRUCTURE without any prior information and then used prior information about different regions (temperate or subtropical) using the LOCPRIOR model 52 . All other STRUCTURE run settings were the same as described above.
Genetic connectivity. To examine the degree of genetic connectivity between temperate and subtropical populations, population genetic structure was examined by AMOVA using Arlequin ver. 3.5 53 . After calculating global F ST , we evaluated the genetic structure using the following prior groupings: (1) all temperate populations vs subtropical populations, (2) recently expanded populations vs pre-existing temperate populations, and (3) populations along the temperate main stream of the Kuroshio Current (ZAH, WAHk, WAHs, KOC, TUK TOI, and KAH) vs populations along the Tsushima Current (AMK, GAHo, and GAHk). Pairwise F ST and pairwise D EST (Jost's direct measure of differentiation 54 ) values across populations were calculated to estimate relatively isolated populations using GenAlEx 6.503 44 . The significance of population differentiation was also determined based on Fisher's exact test using Genepop on the web ver. 4.0.10 45 with Markov chain parameters with 10,000 dememorization steps, 1,000 batches, and 10,000 iterations per batch. Sequential Bonferroni correction was applied to determine statistical significance for all pairs of populations (α < 0.05).
GeneClass2 55 was used to estimate recent gene flow between temperate and subtropical populations using the criteria and the computational algorithm of Rannala and Mountain 56 with default settings. Individuals were assigned to a source population when an individual had the highest assigned probability >70%, and others were regarded as unassigned. Because the two geographically close Wakayama populations (WAHs and WAHk) and the two Goto populations (GAHo and GAHk) were not significantly differentiated (Table 3), we combined these populations for only the GeneClass2 analysis.
IBD. To examine whether neighbouring populations have stronger gene flow than distant populations, the IBD patterns based on pairwise F ST values were regressed onto (1) Euclidean distances between sites, (2) oceanographic distances avoiding physically impossible routes (e.g., passing on land) that were manually estimated using Google Earth and (3) simulated number of accumulated larvae based on a particle tracking model with iterations after 1000 generations (see below, Suppl. 10). We used two different datasets, one including all green lineage populations and the other excluding four recently expanded populations. Significance was tested using the Mantel permutation test (10,000 randomizations) using R ver. 3.4.2. Because the TUK population was surrounded by artificial structures that prevented larval dispersal, this population was excluded from the IBD analysis.
Lagrangian particle tracking simulation. To estimate the probability of larval dispersal between different sites based on ocean current systems, we conducted Lagrangian particle tracking simulation using the Connectivity Modelling System (CMS) 57 with Global HYCOM (Hybrid Coordinate Ocean Modelling 58 ) ocean-current analysis/reanalysis products (https://hycom.org). CMS can evaluate a full range of transport and fate variabilities at high computational efficiency. In the particle tracking model, larvae are treated as passive particles transported on the surface by horizontal current velocities. The particles regarded as larvae were totally released from 111 habitat zones of 1/12° × 1/12° rectangles (~10 km 2 ) (Fig. 6a). The timings of spawning (particle release) from each site were assumed to have begun on the first full moon before the integrated water temperature anomaly from 13 °C from February 1st exceeded 1000 °C day (∫max(T-13, 0)dt). The parameters for estimating the spawning dates (1 °C, February 1 st , and 1000 °C·day) were determined by trial-and-error adjustment using previous records of Acropora spawning timing in Japan. The estimated spawning dates were in general agreement with the records of previous spawning timing (Aizawa et al. unpub. data). The time-series dataset of sea surface temperature used to estimate spawning timing at each site was obtained by extraction from the Global HYCOM analysis/reanalysis products. Spawning events were assumed to take place 3 days before the spawning full moon and to continue for 10 days. From each site, 50 particles were released each day, and a total of 500 particles were released in one year. The maximum larval dispersal period was set for 42 days, which is the maximum length of pelagic larval duration of A. hyacinthus. It was assumed that larvae could begin settling on coral habitats three days after particle release. When the particles passed by the coral habitats, the particles were counted as successful recruitments (the number of larvae dispersed from one place). Finally, a connectivity matrix was obtained by conducting the particle tracking simulation for a total of 24 years from 1993 to 2016. Using the average of the 24-year connectivity matrix, we iterated particle dispersal to simulate the long-term accumulation of dispersed larvae and consider the stepping-stone-like connectivity. This calculation was repeated until the number of the accumulated particle stabilized over time (1000 generations).

Data Availability
The datasets generated during and/or analysed during the current study are available in the Dryad repository (https://doi.org/10.5061/dryad.5ps117g).