Gateways to Hawai‘i: Genetic Population Structure of the Tropical Sea Cucumber Holothuria atra

Holothuria atra is one of the most common and widest ranging tropical, coral reef sea cucumbers in the world, and here we examine population genetic structure based on mitochondrial COI to aid in determining the appropriate scale for coral reef management. Based on SAMOVA, AMOVA and BARRIER analyses, we show that despite its large range, H. atra has hierarchical, fine-scale population structure driven primarily by between-archipelago barriers, but with significant differences between sites within an archipelago as well. Migrate analyses along with haplotype networks and patterns of haplotype diversity suggest that Hawai‘i and Kingman reef are important centers of the genetic diversity in the region rather than an evolutionary dead-end for migrants from the Indo-Pacific. Finally we show that for H. atra Kingman Reef is the most likely stepping stone between Hawai‘i and the rest of the Pacific, not Japan or Johnston Atoll as previously presumed. Based on our data, Johnston Atoll can instead be seen as an outpost of the Northwestern Hawaiian Islands rather than a gateway to the Hawaiian Archipelago.


Introduction
Echinoderms play a major role in structuring many marine ecosystems, and many are described as "keystone species" because of their profound influence on benthic community structure (e.g., [1][2][3], reviewed by Uthicke et al. [4]). In addition to their important ecosystem functions, many echinoderm species are also the focus of artisanal or commercial fishing efforts, particularly the sea urchins and sea cucumbers [5][6][7]. The influence of echinoderm harvest on a wide range of other commercial fisheries, such as abalone, lobster, kelp, and kelp-associated fin fish, has long stimulated discussions of multispecies approaches to managing their exploitation (e.g., [5], reviewed by Purcell [7]). Delineation of the appropriate spatial scales for management zones within a spatial management network requires a detailed understanding of dispersal pathways and population connectivity (reviewed by Hedgecock et al. [8], Thorrold et al. [9], Fogarty and Botsford [10]).
Understanding connectivity in the sea is complicated by the fact that many marine organisms share a biphasic life cycle typified by an adult form that is relatively sedentary and a larval form that can potentially disperse across large expanses of open ocean [11][12][13][14][15]. For example, in the sea urchin genus Tripneustes, some well-known biogeographic barriers, such the Isthmus of Panama or the long stretch of deep water in the western Atlantic, are important barriers to dispersal whereas others, such as the Eastern Pacific Barrier, show no evidence for limiting dispersal [16]. However, the geographic limits of such dispersal are uncertain because it is virtually impossible with current technology to directly track these microscopic juveniles during the pelagic phase (reviewed by Levin [17]) making indirect methods of quantifying larval dispersal particularly attractive (reviewed by Hedgecock et al. [8], Grosberg and Cunningham [18], Selkoe et al. [19], and Hellberg [20]). Proxies for dispersal, such as pelagic larval duration (PLD) and geographic range, have generally been used as rules of thumb in the absence of a detailed understanding of connectivity for most marine species. Unfortunately, intuitive expectations of larval dispersal potential as a function of PLD and 2 Journal of Marine Biology range size are not upheld in recent meta-analyses of the existing literature [21][22][23][24][25][26]88]. Realized dispersal distance is typically less than potential dispersal distance because of the presence of biophysical or biogeographical barriers [27][28][29][30][31]. Barriers that limit dispersal between marine populations include obvious geographical features such as land masses like the Isthmus of Panama [32], but also more subtle factors such as currents and oceanographic regimes [33][34][35][36][37]. The correlation between geographic distance and the probability of larval exchange among sites is low in many marine systems (e.g., [38]), and thus quantitative estimates of connectivity are an important prerequisite for delineating the appropriate scale over which marine populations ought to be managed.
The Hawaiian Archipelago lies at the periphery of the tropical Central Pacific and is the most isolated island chain in the world, making it biogeographically partitioned from the rest of the Pacific Islands (reviewed by Ziegler [39]). This isolation results in one of the highest proportions of endemism in the world (e.g., [40][41][42]; reviewed by Ziegler [39], Eldredge and Evenhuis [43]). Though there are many examples of pan-pacific coral reef organisms in Hawai'i, the isolation of the Hawaiian Archipelago is thought to limit larval exchange sufficiently that colonization is rare [44]. For example, Kay [45] estimated that Western Pacific marine species successfully colonize the Hawaiian Archipelago about once every 13,000 years. Unlike the terrestrial fauna, however, the Hawaiian marine fauna contains a large proportion of endemics that are differentiated but not diversified from their Indo-West Pacific roots [39,[46][47][48]. Johnston Atoll is believed to be a stepping stone into Hawai'i, and simulations of larval dispersal suggest that larvae from Johnston atoll can reach French Frigate Shoals or Kaua'i along two separate larval corridors [49,50].
The lollyfish, Holothuria atra, is one of the most common shallow-water tropical sea cucumbers in the Indo-Pacific, spanning from Madagascar to French Polynesia [51,52]. H. atra performs vital ecosystem services on coral reefs for which there is an active fishery in many regions of the Pacific [7,53,54]. Echinoderms are described as a boom-bust phylum in which populations go through marked natural population cycles [4], an attribute that can compound problems in a harvested population but may hasten repopulation in previously impacted areas. As such, there is a call for ecosystem-based management of sea cucumber harvests [7]. Furthermore, the boom-bust nature of echinoderms has important implications for connectivity in evolutionary time frames, where biological attributes can drive population structure to a greater extent than oceanographic processes as hypothesized in the Tripneustes sea urchins [16]. Together these characteristics make H. atra an ideal organism to examine levels of connectivity and historical population dynamics to inform management and to test hypotheses about Hawai'i's connection with other archipelagos in the Central Pacific. Here, we assess the inferred range of dispersal for H. atra in Hawai'i and the Central Pacific by investigating its mitochondrial genetic population structure in an attempt to delineate the appropriate scales for management.

Methods
2.1. Sampling, PCR, and Sequencing. Holothuria atra was sampled from five archipelagos (Hawaiian Islands, Line Islands, Marshal Islands, Bonin Islands, and Ryukyu Islands) at a total of 19 sites (Figure 1). Sampling in the Northwest Hawaiian Islands and the Line Islands took place on research cruises aboard the NOAA R.V. Hi'ialakai. All other samples were collected on shore dives or while snorkeling. Sampling took place between spring 2006 and fall 2009. Samples were obtained nonlethally through muscle-tissue biopsy and preserved in either 95% ethanol or DMSO salt buffer and archived at the Hawai'i Institute of Marine Biology at room temperature. Skillings and Toonen [55] made an extended discussion of sampling and preservation protocol. No asexual morphs-distinguished by transverse scarring, smaller body size, and their location in lagoonal habitats-were found during sampling expeditions, and no reports are known indicating the presence of the asexual stage of H. atra in the sampled locations. The asexual morph of H. atra appears to be located only in the Southern and West Pacific (e.g., [52,56,57]).
Total genomic DNA was extracted using DNeasy Blood and Tissue Kits (QIAGEN) following the manufacturer's instructions. Polymerase chain reaction (PCR) was used to amplify a 423-base pair fragment of the mitochondrial cytochrome c oxidase subunit I gene (COI) using custom primers created with Primer3 [58] targeting Holothuria spp.: GenHol2L (5 -AACCAAATGGTTCTTGCTTACC-3 ) and GenHol2R (5 -TTCTGATTAATCCCACCATCC-3 ). PCR was performed using 15 μL reactions containing 1 μL of diluted DNA extract (one part template DNA to 199 parts nanopure water), 1 μL each of 0.2 μM forward and reverse primers, 0.6 μL of 0.5 μM BSA, 7.5 μL of (Bioline) Biomix Red diluted as per manufacturer's instructions, and 3.9 μL of nanopure water. PCR was done on Icycler thermocyclers (Bio-Rad Laboratories) with an initial denaturation at 95 • C for 7 min followed by 35 cycles of a denaturing step at 95 • C for 1 min, annealing at 50 • C for 1 min, and extension at 72 • C for 1 min. A final extension at 72 • C was held for 7 min before refrigeration. PCR product (8 μL) was treated with 0.7 μL of Exonuclease I combined with 0.7 μL of calf intestinal alkaline phosphatase (Exo-CIAP) and incubated at 37 • C for 30 minutes and with a final inactivation step at 85 • C for 10 minutes. The treated PCR product was sequenced using an ABI Prism automatic sequencer at the Hawai'i Institute of Marine Biology's EPSCoR sequencing facility. All samples were sequenced in the forward direction; uncertain sequences and all unique haplotypes were also sequenced in the reverse direction for confirmation. Sequences were compiled and trimmed using Sequencher 4.8 and aligned using ClustalW implemented in Bioedit 7.0.5 [59,60].

Data Analysis.
A statistical parsimony network of mitochondrial haplotypes was constructed by creating a reduced median network that was then used to make a median joining network; both procedures implemented in Nei's average pairwise genetic difference (π) [63] and haplotype diversity (h) were calculated in DnaSP 4.1 [64]. The effective number of alleles was calculated by hand following Jost [65]. Tajima's D [66] and Fu's F S [67] were calculated in Arlequin v.3.1 with 10,000 simulations to establish confidence intervals.
To assess levels of genetic differentiation between sites, we calculated pairwise Φ ST values using Arlequin 3.1 [68] and pairwise D est chao values by hand following Jost [65]. Φ ST is a fixation index incorporating genetic distance that ranges from 0 to 1, where a zero indicates identical haplotypic composition and a one signifies alternate fixation of alleles and a complete lack of gene flow. D est chao is an index of genetic differentiation, which does not account for genetic distance among haplotypes but also ranges from 0 to 1 (note that both Φ ST and D est chao can be slightly negative due to bias correction for sampling error). In the case of D est chao , a zero also indicates identical haplotypic composition, but unlike Φ ST , a one simply indicates that no haplotypes are shared between the populations. The primary difference in interpretation is that in the absence of gene flow Φ ST values can be significantly less than one, while this is not the case for D est chao , which is argued to be an advantage of this latter statistic [65]. To correct the critical P value for statistical significance in pairwise comparisons, the familywise false discovery rate (FDR) correction found in [69] was implemented. Analysis of molecular variance (AMOVA) was used for hierarchical analysis of the partitioning of COI diversity among sites within archipelagic regions and among archipelagic regions using Arlequin 3.1. SAMOVA 1.0 was used to identify groups of samples that maximize the proportion of total genetic variance due to differences between regions [70]. The most important genetic barriers were ranked using BARRIER 2.2 [71]. BARRIER uses Monmonier's maximum-difference algorithm to compare a matrix of difference values, such as pairwise Φ ST values, with a matrix of geographic distances in order to identify the strongest barriers within the matrix. We compared barriers created using each Φ ST and D est chao distance matrix. AMOVAs were performed using groupings determined by SAMOVA and BARRIER for hypothesis testing to compare genetic groupings to the archipelagic groupings. The pairwise Φ ST and AMOVA analyses were conducted using a distance matrix with 50,000 permutations and the Tamura-Nei mutational model [72] with gamma = 0.0164. The mutational model HKY+G was selected using AIC in MODELTEST 3.7; the model hierarchy was used to select the closest available model when the best-fit model could not be implemented by the chosen program, as in the case of ARLEQUIN [73]. Regardless, the inferences are robust to the mutational model, and our conclusions are not altered regardless of which model is chosen (data not shown). MRBAYES 3.1 was used to construct a Bayesian estimation of a phylogeny containing all H. atra haplotypes from this study along with all Holothuria and Actinopyga COI haplotypes available from GenBank as of February 2010 [74]. Two independent runs with identical conditions were completed and averaged. A general time reversible (GTR) simple nucleotide model with a gamma-shaped rate variation of 0.0164 was used; Markov chain length = 4 × 3,000,000 sampled every 100 generations with a 10% burnin. The GTR nucleotide model was chosen as it is the most general and neutral nucleotide model available in MrBayes 3.1 and corresponds most closely to the Tamura-Nei model [75]. The sea cucumber Actinopyga agassizi was set as the outgroup. Program defaults were used for all other settings. MrBayes was used to summarize all of the trees produced into a single consensus tree.
RAxML 7.0 [76] implemented through CIPRES Web Portal v.1.15 [77] was used to construct the highest scoring maximum likelihood-based estimation of a phylogeny containing all haplotypes used in the Bayesian analysis and run 10,000 bootstrap simulations to assess branch support. A GTR nucleotide model that uses four discrete gamma rates set by the program was used for the analysis; program defaults were used for all other settings.
Bayesian coalescent-based calculations of migration rate among regions (N e M) and the region mutation parameter (θ) were conducted using MIGRATE 3.1.3 [78]. Three independent runs of a Bayesian MCMC search strategy were completed and averaged by MIGRATE. A nucleotide model with a transition-to-transversion ratio of 6.1584 : 1 and three regions of substitution rates with a gamma-shaped rate variation of 0.016 was used; Markov chain length = 1,000,000 sampled every 20 generations with a 10% burn-in. Program defaults were used for all other settings. The transition-totransversion ratio was calculated using Modeltest 3.7. Two replicate MIGRATE analyses were run using different population groupings. Preliminary analyses that split the data by sampling location returned flat posterior probabilities, presumably from having too many parameters to estimate. The software's author advocates using the minimal number of sensible regions in order to reach convergence (Peter Beerli, pers. comm.). The first analysis used regions separated along the most important breaks identified by the program BARRIER. Archipelagos were used as regions for the second analysis with the Hawaiian Archipelago further divided between the main Hawaiian Islands and the Northwestern Hawaiian Islands plus Johnston Atoll and Kingman Reef separated from the Line Islands; this grouping was chosen based on BARRIER and pairwise analyses. Values for the migration rate among regions (M) and region mutation parameters (θ) were taken from the highest peaks in the posterior probability distribution curves. The posterior probability distributions were examined to determine the validity of each estimated parameter.

Results
A total of 385 individuals, 55 haplotypes, and 37 private haplotypes are sampled in this study (Table 1). Of the 18 haplotypes shared across sampling sites, 10 are found in multiple archipelagos ( Figure 2). However, no haplotypes are shared between the most distant regions: Hawai'i and Japan. Because many population genetic estimates are relatively insensitive to weak selection [79], loci which do not show significant deviations from neutral expectations should provide reliable inferences about population structure [80]. None of the site-by-site Tajima's D values were significant, and only Laysan deviated from expectation using Fu's F s ; thus, there is no evidence to indicate that nonneutral processes are responsible for the pattern of COI haplotype diversity presented here.
To ensure there was no misidentification of the samples included here, we performed a phylogenetic reconstruction of our samples with those available in GenBank. We confirmed that the samples included here are monophyletic and there are no reciprocally monophyletic groups among the H. atra haplotypes included in our population genetic analyses (Appendix A).
Haplotype diversity, as a function of longitude, increases from west to east across the Pacific Ocean (Table 1, R 2 = 0.80, P < 0.05). Japanese sites exhibit the lowest haplotype diversity (h = 0.51-0.66) and effective number of haplotypes (H E = 2.0-2.9); the centrally located Line Islands (h = 0.65-0.90, H E = 2.9-9.7) and Marshall Islands (h = 0.75-0.76, H E = 4.0-4.1) exhibit mid to high levels of diversity; the Hawaiian sites exhibit the highest diversity (h=0.75-0.90, H E = 4-10). Nucleotide diversity does not appear to be correlated with haplotype diversity because nucleotide diversity is lower in the Hawaiian Archipelago than in all other locations except for Okinawa (Table 1). This pattern can be visualized in the haplotype network, where sites exhibiting high nucleotide diversity harbor disparate haplotypes separated by a relatively large number of mutations ( Figure 2).
Four AMOVAs were run on the H. atra COI haplotype data ( Table 2). In each AMOVA, a different method was employed to group the population samples into regions. Under the first geographic hypothesis, sites were grouped into five regions by archipelago, with Johnston Atoll included with Hawai'i. In order to assess the subdivision of the Hawaiian Archipelago into the NWHI and the MHI, this grouping hypothesis was compared to a second geographic hypothesis, where population samples were grouped into six regions; Hawaiian sites were divided into the MHI and the NWHI + Johnston and four regions comprised of the four remaining archipelagos. Both grouping hypotheses (geographic hypothesis one and two) exhibited similar partitioning of variation among groups (32.1% versus 30.9%) and among populations within groups (8.9% versus 5.5%). In both cases, there was stronger partitioning among the groups of samples (Φ CT = 0.32, Φ CT = 0.31; P < 0.0001) than among the simples nested within the groupings (Φ SC = 0.13, Φ SC = 0.07; P < 0.0001).
The two geographic hypotheses were compared to sixregion groupings identified by SAMOVA and BARRIER; BARRIER selected the same dominant barriers using both the Φ ST and D est chao distance matrices (Table 2 and Figure 1). Six regions were chosen for a direct comparison to the archipelagic geographic hypothesis that included Hawai'i divided into two main regions. These groupings partitioned variance similarly to the geographic hypotheses with the SAMOVA grouping minimizing among-population withingroup variance with more among-group variance explained ( Table 2). The grouping of population samples using BAR-RIER and SAMOVA had slightly greater levels of genetic Figure 2: Haplotype network for Holothuria atra. Each circle represents a unique haplotype connected by a line to those that differ by one base pair. Nodes on lines indicate a missing haplotype, and numbers represent multiple missing haplotypes. Each haplotype is color coded by site, and circle size is proportional to frequency. The smallest circles represent one occurrence of a halotype.  The results from the MIGRATE runs show similar patterns of gene flow between regions (Table 4, Figure 1). Effective migration rates (N e M) between regions are low. There is less than one migrant per generation, the ruleof-thumb number below which population cohesion starts to break down, between most regions [81]. The exceptions include the one-way migration from the Main Hawaiian Islands into the Northwest Hawaiian Islands and Johnston Atoll and the one-way migration from Kingman Reef to the Main Hawaiian Islands ( Table 4). The high effective migration rate from the NWHI to the MHI in the first analysis splits almost evenly between Kingman Reef and the MHI when these two regions are separated in the second analysis (Table 4). Overall, higher effective migration rates are observed leaving the Hawaiian regions then going into them (Table 4, Figure 1). Though effective migration rates are a product of migration and effective population size, the effective migration rates larger than 1 migrant per generation are driven primarily by migration and not effective population size (Appendices A and B). This pattern is indicative of recent migration rather than ancestral polymorphisms and high effective population sizes. Posterior probability distributions for all values were in the form of unimodal curves. A full description of N e and M values for both analyses can be found in Appendices B and C.

Discussion
In this survey of population genetic structure, we elucidate patterns of connectivity throughout the north-central range of the sea cucumber Holothuria atra with a focus on the Hawaiian Archipelago. The Hawaiian Archipelago is highly isolated and also contains one of the highest proportions of endemism in the world (e.g., [40][41][42]; reviewed by Ziegler [39], Eldredge and Evenhuis [43]). Though there are many pan-pacific marine organisms in Hawai'i, the isolation of the archipelago is thought to limit larval exchange such that colonization is rare but sufficient to maintain species cohesion among these taxa. The Hawaiian marine fauna contains a large proportion of endemics that are differentiated but not diversified from its Indo-West Pacific roots [39,[46][47][48]. In this scenario, Hawai'i is seen primarily as a dead end, an isolated land mass that does not contribute in a significant way to the overall diversity of the tropical pacific.
Counter to the island biogeography hypotheses of Hawaiian diversity, Jokiel and Martinelli [44] proposed the Vortex model of speciation, wherein the stunning biodiversity of the Coral Triangle is a result of centrifugal accumulation of species from the peripheral habitats around the Pacific. Though these two models primarily make predictions about speciation-level processes and do not speak directly to gene flow within a species, they do make opposite claims about the dominant direction of gene flow and dispersal. H. atra has a broad species range, extending from the Western Indian Ocean to the Eastern Pacific Ocean, which suggests the capacity for long-distance dispersal; however, populations showed significant population structuring within archipelagos, sometimes across very short oceanic distances. Even so, hierarchical genetic population structure in H. atra gives insight into the phylogeography of the north-central tropical Pacific. Our data test between the divergent hypotheses of whether peripheral archipelagos act as a source of genetic diversity in the Pacific and the likely colonization routes, into and out of, the extremely isolated Hawaiian Archipelago.

Biogeography and Range
Size. If a large species range is a consequence of high dispersal potential, then H. atra should have little pronounced population structure, especially across small scales [11,15,82]. Indeed, this is the case for many species in the central West Pacific [16,[83][84][85][86]. Despite a species range which stretches from the Western Red Sea to the eastern Central Pacific in which H. atra is found in almost all shallow tropical habitats, we did not find support for extensive dispersal. The majority of sites from which we sampled H. atra were genetically distinct, with some sites less than 75 km apart being among the most distinct in our study (Table 3). These contrasting patterns highlight the dangers of making predictions about population connectivity and diversity based solely on the location and size of a species' range.
The larval life history of H. atra is not known exactly, but they require at least 18-25 days to reach competency to settle and are capable of traversing long oceanic distances with sufficient frequency to maintain species cohesion across a very broad geographic range [87]. The obvious question becomes why then is population subdivision found on such small geographic scales (e.g., Kingman Reef and Palmyra Atoll are only 67 km apart)? Counter to intuition, the geographic distance among sites is a poor predictor of the ease with which larvae can disperse among locations; the "oceanographic distance" experienced by larvae between sites is uncorrelated with geographic separation between them [36,38]. Likewise, recent meta-analyses indicate the relationship between the length of pelagic larval development and dispersal ability is not as tight as has been generally assumed [23][24][25][26]88]. Finally, a broad meta-analysis by Lester et al. [21] indicates that the intuitive relationship between range size and larval dispersal potential is poorly correlated overall but can play an important role in some taxa. Toonen et al. in this issue also show a number of breaks in the Hawaiian Archipelago that are shared by several species and are unexplained solely by appeal to one metric such as range size or larval dispersal potential. Although the mechanism of isolation across small scales remains unknown, our data clearly indicate that H. atra is not one of those species for which range size predicts relative dispersal ability.

Population Structure in the Hawaiian Archipelago and
Johnston Atoll. Our mtDNA examination of Holothuria atra reveals significant genetic population structure across the surveyed portion of the range. There are two interesting patterns to this structure. Excluding Laysan Island, there are no significant pairwise differences between any other islands in the NWHI (spanning nearly 2000 km), suggesting that the NWHI, excluding Laysan, comprises a single large population. In contrast, there is significant structuring within the MHI (roughly 600 km) and between the NWHI and the MHI. This finding suggests that factors beyond merely geographic distance influence population partitioning.
Johnston Atoll, the nearest neighboring land mass, roughly 860 km south of French Frigate Shoals, is genetically distinct from most of the MHI and Laysan and genetically similar to all of the NWHI except Laysan. It has been suggested that Johnston Atoll acts as a stepping stone into the Hawaiian Islands [89]. Kobayashi [49,50] used computer simulations to predict two larval transport corridors from Johnston Atoll to the Hawaiian Archipelago: one corridor stretching from Johnston to French Frigate Shoals in the NWHI and one from Johnston to O'ahu in the MHI. Our data support the predicted larval transport corridor between Johnston Atoll and French Frigate Shoals, but not the corridor predicted between Kaua'i and Johnston. Additionally, based on our data, Kingman Reef may also be an important stepping stone into and out of Hawai'i. The BARRIER analysis shows the division between the NWHIs, including Johnston Atoll, and the MHI to be the strongest barrier to gene flow within the Archipelago (Figure 1). Migration across this barrier is heavily one sided, where migration from the MHI into the NWHI dominates. The effectively one-way migration rates into the NWHI and Johnston Atoll coupled with the strong genetic similarity between Johnston Atoll and the NWHI suggest Johnston Atoll is an isolated outpost of the Northwest Hawaiian Islands, providing support for a vortex model [44] rather than the stepping stone entry into Hawai'i [89] for H. atra.
These data indicate that Johnston Atoll exchanges migrants with Hawai'i far more often than its nearest neighbors to the south, and the same can be said for Kingman Reef. This result is particularly surprising because in the case of Kingman Reef, H. atra sampled there show greater similarity to populations in Hawai'i (roughly 1700 km southwest of Honolulu) than they do to those sampled at Palmyra Atoll, only 67 km away.

Phylogeographic Relationships between Archipelagos.
Counter to conventional wisdom that Hawai'i is a passive recipient of rare dispersal from the diverse Pacific, the weight of available evidence, including pairwise Φ ST values, mtDNA phylogeny, BARRIER divisions, and clustering within the haplotype network, provides substantial evidence for the opposite pattern in H. atra; Johnston Atoll is an outpost of Hawaiian diversity, and Kingman Reef acts as the primary stepping stone between the Hawaiian Archipelago and the rest of the Pacific. As far as we are aware, this is the first time empirical evidence has been provided for such a pathway. Higher haplotypic diversity in Hawai'i and the Line Islands relative to the other archipelagos supports a scenario in which population sizes are far greater, or Hawai'i and/or Kingman Reef are the ancestral population in the region. Likewise, the dominant haplotypes found in the Japanese Archipelagos are relatively distantly related and appear derived (Figure 2), suggesting that the western portion of the surveyed range was colonized in at least two separate events (or one of them has gone extinct in Hawai'i and Kingman), one of which did not make it all the way to Okinawa.
Excluding the dense sampling within the Hawaiian Archipelago, virtually all pairwise comparisons between sites are significantly different from each other ( Table 3). The few comparisons that were not significant between sites have relatively high pairwise values but low sample sizes, a likely statistical limitation also noted by Bird et al. [90]. Hierarchical population structuring was detected with AMOVA using either D est or Φ ST values; sampling sites within archipelagos are significantly different from each other but are more similar within than between archipelagos. Four of the five most substantial restrictions to gene flow uncovered in this study (and the top ranked by BARRIER) were those between archipelagos; the one exception to this trend being Kingman Reef which is included with the Hawaiian rather than the Line Islands (Figure 1). The AMOVA run using the regions selected by BARRIER was only minimally different than the AMOVA run using regions divided by archipelagos; grouping Kingman Reef with the MHI, as per BARRIER, did explain 1.5% more of the overall variance.
In addition to the distinct archipelagic groupings in the haplotype network (Figure 2), several other patterns are noteworthy here. First, the NWHI and Johnston Atoll haplotypes are clustered together and interspersed whereas the MHI haplotypes are clustered together. Also, the Japanese haplotypes occur in two divergent areas of the network. The Line Island and Marshall Island haplotypes are interdispersed throughout the network, suggesting that these island groups are either mixing or transition zones. Nearly every locality haplotype (those found in only one sampling location) branch off in a starburst pattern from the major haplotypes found primarily in the same archipelago. This pattern is an indication that regional populations have been separated long enough for new haplotypes to arise, and that these new haplotypes are not being spread to other archipelagos by long-distance dispersal. Uniformly low migration rates between archipelagos estimated with MIGRATE support this isolation scenario.

Conclusion
Many echinoderm species are the focus of artisanal or commercial fishing efforts, and managing these fisheries requires a detailed understanding of dispersal pathways and population connectivity within a spatial management network. The Hawaiian Archipelago lies at the periphery of the tropical Central Pacific and is the most isolated island chain in the world; the question remains as to why some species maintain connectivity and species cohesion between the Hawaiian Islands and the rest of the Pacific, why some species diverge and become Hawaiian endemics, and why other species with similar inferred dispersal ability fail to colonize the Hawaiian Archipelago at all.
The genetic diversity of COI in H. atra across the studied portion of the range presents a complex pattern, but it is not inscrutable. Based on AMOVA, SAMOVA, and BARRIER analyses, it can be seen that population structuring is hierarchical; there are significant differences between sites, but the primary degree of population structure is archipelago by archipelago. Our analyses taken together suggest that the Hawaiian Archipelago and Kingman reef are ancestral populations in the region with migration moving out of these periphery archipelagos toward a less diverse central Pacific rather than the reverse. This pattern is inconsistent with the hypothesis that Hawai'i is a dead end for rare migrants from the Indo-Pacific. Instead, the weight of the evidence shows that these peripheral populations are not sinks, but important centers for the generation of genetic diversity feeding back towards the West Pacific. Specifically for H. atra, our data suggest that the pathway between Hawai'i and the rest of the Pacific is primarily out through Kingman Reef and the Line Islands and not in through Japan, the Marshall Islands, or the closest neighbor to the Hawaiian Archipelago, Johnston Atoll. We show that, at least for H. atra, Johnston Atoll is in fact an outpost of the Northwestern Hawaiian Islands and not a primary gateway for colonization of the Archipelago.
Considerable evidence is accumulating that it is indefensible to make predictions of connectivity based solely on proxies such as ecological or phylogenetic similarity, pelagic larval duration, or species range sizes [21,[23][24][25]90]. The fine-scale structuring of populations in H. atra suggests that place-based management approaches, as exemplified by ecosystem-based management, are ideal for responding to the complex relationships between genetically distinct populations. Holothuria atra must be managed on a local scale; migration between archipelagos, and often between islands, does not occur in ecologically relevant time frames.

A.
Phylogenetic tree of sampled Holothuria atra haplotypes and Holothuria and Actinopyga haplotypes retrieved from GenBank. The maximum likelihood tree is shown; Bayesian analysis produced a tree with almost identical topology. The first of the numbers beside branches corresponds to bootstrap support calculated from the maximum likelihood analysis; the second number, seperated by a /, corresponds to the posterior probability estimated during the Bayesian analysis. Only support values with 50% majority rule or 0.50 posterior probability are included. Actinopyga agassizi is the outgroup. All non-H. atra nodes were collapsed. GenBank accession numbers are included for uncollapsed taxa (see Figure 3).