Population structure among octocoral adults and recruits identifies scale dependent patterns of population isolation in The Bahamas

Sample collection

In order to assess genetic diversity and structure of Antillogorgia elisabethae tissue samples from adult colonies were collected from 13 sites in the northern Bahamas (Table 1 and Fig. 1). Colony samples were collected from arbitrarily chosen colonies at each site. There have not been any observations of reproduction via fragmentation in detailed surveys for recruits at multiple sites in The Bahamas (Lasker, 2013). Colonies were collected at several meter or greater intervals in order to minimize the likelihood of sampling familial groups. Collections at some sites required searching along paths of 200–400 m. Depth at the different sites ranged from 8 to 22 m. From each colony, ca. 4 cm of live tissue was removed and preserved in 95% ethanol. Many of the sites were similar to those reported in Gutiérrez-Rodríguez & Lasker (2004c), but only two, (San Salvador and Hog Cay) were identical sites and samples.

Recruits were collected from a series of 4 sites (Table 1) along the Little Bahama Bank (LBB) every year between 2004 and 2012, with the exception of 2008. At each site all of the recruits within 20 randomly placed 1 m² quadrats were collected and preserved in 95% ethanol. Details of the sites and quadrat placement are presented in Lasker (2013). Collections prior to 2007 were made in May and early November, respectively 5 and 11 months after November/December spawning. In order to collect recruits closer to the time of spawning (and thus subject to less post-settlement mortality, the sampling was changed to January and May/June in 2009 and thereafter. Recruits collected in January were usually single polyps. In both sampling schemes the vast majority of recruits were collected in the first of the two months. All collections were conducted under a research permit from the Department of Marine Resources, The Bahamas.

Laboratory analysis

Allelic variation was assessed by genotyping a total of 704 samples of A. elisabethae adults and 776 recruits. Genomic DNA was isolated from adults or large recruits using 0.5 cm of tissue and from the entire sample among smaller recruits. The CTAB phenol-chloroform-isoamyl alcohol protocol (Coffroth et al., 1992) or DNeasy Tissue Extraction kit (Qiagen) was used for extractions. Loci used in the analyses and references for their use were Pel19, Pel32 (Gutiérrez-Rodríguez & Lasker, 2004a; Gutiérrez-Rodríguez & Lasker, 2004b), Pel62, Pel84 (Lasker et al., 2008) and Ael-AG-1, Ael-AGC-1, Ael-ATC-1, Ael-AGAT-1, Ael-AGAT-2, and Ael-AGAT-3 (Porto-Hannes & Lasker, 2013). In brief, PCR consisted of 10-µL reactions, containing 2 U of Taqpolymerase, 10 mM Tris–HCl (pH 8.3), 50 mM KCl, 200 µM dNTPs, approximately 10 ng of DNA template, 1.5–2.5 mM MgCl₂, and 0.025–0.05 µM of each primer. Cycling conditions were 2 min at 94 °C, followed by 28–32 cycles at 94 °C for 30 s, 50–56 °C for 45 s and 72 °C for 45 s, with a final extension period of 10 min at 72 °C.

PCR products were screened on 7% polyacrylamide gels in a LI-COR NEN^® Global IR2 DNA Sequencer System, using fluorescently labeled primers as described in Santos et al. (2003). We determined allele sizes by comparing the migration of the PCR products to 50-350 bp size standards (LI-COR Biotechnology Division) regularly spaced across the gel. Scoring was manual with allele sizes frequently verified by comparing migration distances of alleles with those of size standards using ImageJ (NIH).

Population genetics and statistical analyses

Adults

A total of 704 adult colonies were analyzed from the 13 sites. Summary data for the loci and sites are presented in Table 2 and Table S1. Numbers of alleles ranged from 1 to 43 at individual sites (Table S1) and across all sites ranged from 9 to 50. Heterozygote deficiency was observed at many loci and populations and global tests of Hardy-Weinberg equilibrium (FSTAT v.2.9.3) identified significant deviations from Hardy-Weinberg equilibrium across both sites (not shown) and loci (Table S2). Microchecker identified 2 loci, Pel-62 and Ael-AGAT-1, which had null alleles at 10 of the sites and one locus, Ael-AGAT-2, having null alleles at 9 of the sites. Unless indicated otherwise those 3 loci were not included in subsequent analyses. After dropping the 3 loci with putative null alleles, 23 of 89 site by locus tests for deviation from Hardy-Weinberg equilibrium were significant (P ≤ 0.001), but there was no consistent pattern across either loci or sites (Table S2). Examining all combinations of loci and populations, there were 38 cases in which loci were not in genetic equilibrium (P ≤ 0.05, Table S3), but none of those cases were significant after probabilities were corrected for multiple comparisons.

Global F_ST calculated in FSTAT v.2.9.3 was 0.177 and that value as well as the vast majority of pairwise F_ST values were significantly greater than 0.0 (Table 3). F_ST calculated with FSTAT v.2.9.3 and FreeNA (Table S4) were very similar to each other. F_ST values indicate the populations from HC (Exuma Sound) and PR (San Salvador) were clearly differentiated from each other and from all of the other sites. The TOTO and AN sites were also dissimilar to most other sites but to a lesser extent. Although both sites were situated along the western edge of the Tongue of the Ocean, AN was no more similar to TOTO than to the other sites. TR, although separated from the other sites by Andros Is. and the Great Bahama Bank, also had intermediate levels of genetic differentiation to both LBB and GBB sites (F_ST = 0.0404–0.1469, Table 3). With the exception of the WC population at the western edge of Grand Bahama, comparisons among the populations on the LBB had low (<0.05) pairwise F_ST values. There were few cases of non-significant F_ST values and 6 of the 7 such cases were among LBB sites and most commonly between adjacent sites. The seventh case was between EGG and GP, which are not adjacent sites and GP is on the LBB, which is separated from EGG by 60 km at its closest. EGG and all of the LBB sites (excluding WB) had low F_ST values, with the adjacent CHSL and CHR having particularly low, but significantly greater than 0.0, values of 0.0280 and 0.0296 respectively. As expected from the nature of the index, $G_{\mathrm{ST}}^{″}$ values were greater than F_ST and are shown in Table S5. One pairwise comparison (GP vs. CHR) in which the F_ST was not significantly greater than 0.0 had a $G_{\mathrm{ST}}^{″}$ that was significantly greater than 0.

The pattern of increasing genetic differentiation with geographic distance suggested by the F_ST values is evident in the Mantel tests (Fig. 2A). Plots of genetic distance versus geographic distance were significant (Fig. 2A, P = 0.01) and the results were very similar regardless of whether geographic distance or adjusted geographic distance was used (geographic distance, y = 0.0007x − 0.0208, R² = 0.6775, P = 0.010; adjusted geographic distance, y = 0.0005x + 0.0008, R² = 0.4857, P = 0.010). The proportion of variance explained by the relationship was greater using geographic distance and those data are presented. As evident in Fig. 2B pairwise comparisons among sites from the LBB as well as those comparing LBB sites with EGG form a clear linear pattern. If the regression for the LBB is recalculated including the comparisons with EGG, the slope of the line is unchanged and the regression explains 84% of the variation. The regression that included comparisons among the other sites (Fig. 2A) and of those sites to LBB sites also generated a significant regression, but inspection of the data suggests there are several subsets of comparisons, that do not appear to change with distance. The pairwise comparisons with F_ST/(1 − F_ST) values >0.15 are all comparisons of Exuma Sound (HC) and San Salvador (PR) with the other sites, and each of those sets of comparisons shows little change with distance and perhaps even an inverse relationship.

The Structure and Geneland analyses of the populations clearly characterized the complex nature of population structure of A. elisabethae in The Bahamas. An initial Structure analysis, which included all 13 sites exhibited a best fit model of K = 2 (data not shown) in which colonies from Exuma Sound (HC) and San Salvador (PR) generally belonged to one population and those from the other sites to a second population. Subsequent analysis without the Exuma Sound (HC) and San Salvador (PR) collections generated highly variable results. LnP[k] values from the Structure analyses did not monotonically increase with increasing K, and there were high levels of variance at K values >3, both of which suggested the simulations were not converging. Increased steps (200,000 burn in and 10⁶ MCMC replications) had no appreciable effect on the results. Repeating the analysis using all 10 loci; i.e., including those in which null alleles were present, also indicated K = 2 when all 13 sites were included (Figs. 3A and 3B). The individuals from San Salvador (PR) belonged to a distinct population with those from Exuma Sound (HC) exhibiting admixture between the San Salvador (PR) population and a second population, which was represented by the colonies from the remaining sites. That analysis suggested K = 3 among the 11 sites (Figs. 3C and 3D).

Structure is not well suited for populations that fit an IBD model (Pritchard, Stephens & Donnelly, 2000), and thus the meaning of the specific value of K = 3 is unclear. However, the effects of IBD are evident in the manner in which admixture changes between sites in the Structure analysis. The sites in Fig. 3D are arrayed from west to east along the southern edge of LBB and then east to west along the GBB. There was a distinct west to east cline in the contribution of the different “source” populations on the LBB, and individuals from the GBB sites differed from those of the LBB, with the population near Bimini (TR) being the most distinct.

Geneland incorporates the geographic location of the samples in its analysis, and the analysis partitioned the data into K = 8 clusters (Fig. 4). The grouping within the LBB sites and EGG changed with different λ values. When λ = 2, all LBB sites except WC formed one group. At λ = 5–100, WC formed a cluster, BS and SC another cluster and the rest of the LBB populations along with EGG formed a third cluster. At λ = 100, the most likely number of cluster was 10, however, two of these clusters were “unsampled” populations. The other 8 clusters were the same as described above. Although the Geneland analysis partitioned the sites into a greater number of populations, the results are concordant with the patterns seen in the Structure plots (Fig. 3). In both analyses there is a core of sites on the LBB. In the Structure analysis, sites on the LBB differed from adjacent sites by relatively small differences in the levels of admixture, while in some Geneland analyses those were further divided into the sites near Great Abaco (CHR, CHSL, GP, MI and sometimes EGG) and those closer to Grand Bahama (SC, BS). WC and TR which had different patterns of admixture in the Structure analysis were also distinct populations in the Geneland analysis.

Principal coordinate analysis (PCoA) of genetic distances (Fig. 5) identified a pattern similar to that hypothesized in Geneland and Structure. San Salvador (PR) and Exuma Sound (HC) and to a lesser extent the Tongue of the Ocean (TOTO) populations were distinct from each other and from the other sites.

Isolation by distance

IBD has been found in a large array of reef species including reef fishes (Buston et al., 2012; Nanninga et al., 2014; Pinsky et al., 2012; Puebla, Bermingham & Guichard, 2009), gastropods (Crandall, Treml & Barber, 2012), octocorals (Andras, Kirk & Harvell, 2011; Ledoux et al., 2010; Yasuda et al., 2014) and scleractinians (Maier et al., 2005; Polato et al., 2010; Porto-Hannes et al., 2015). The pattern of genetic differentiation of A. elisabethae along the LBB provides a remarkably good fit to an IBD model, explaining 86% of the variation in genetic differentiation along the LBB. The southern edge of the LBB does not support a continuous well developed reef community. None-the-less, sites that support A. elisabethae are apparently common enough to create stepping stone populations, which in concert with A. elisabethae larval dispersal generates an IBD model.

A. elisabethae on the LBB are distributed in a linear array of populations. In such cases, the slope of the regression (b = 2 × 10⁻⁴, Fig. 2B and Fig. S1) between geographic distance and F_ST/(1 − F_ST) can be used to estimate variance of the dispersal kernel using the equation b = 1/(4Dσ²), where D is the population density and σ² is the variance of the dispersal kernel (Rousset, 1997). Density of reproductive colonies from 8 sites between Great Abaco and Grand Bahama (CHSL to SC), at which the species was abundant, averages 0.74 m⁻² (H Lasker, 2011, unpublished data). Assuming the suitable habitat zone is 200 m wide the estimate of density is D = 1.48 × 10⁵ km⁻¹ of reef tract. If A. elisabethae was distributed continuously that would correspond to σ of 83 m. However, A. elisabethae is patchily distributed along the length of the LBB, with several kilometers of high abundance often followed by kilometers in which the species is virtually absent. In addition, all reproductive size octocoral colonies are not equally successful in either producing or fertilizing eggs (Beiring & Lasker, 2000; Lasker et al., 2008), which further reduces the effective population size. Using an effective density only 1% of the density observed at the study sites only increases the estimate of σ to 830 m. The relatively small variance of the dispersal kernel suggests a high level of self-recruitment, which is consistent with the assignment based analyses. This conclusion is also concordant with an analysis of recruitment rates at these and 2 additional LBB sites. In that study the density of reproductive colonies within a 100 m diameter area at each of the 8 sites accounted for 46% of the variation in recruitment between the sites, again suggesting very limited dispersal (Lasker, 2013).

In general, brooding corals exhibit greater population genetic structure and small inferred dispersal distances than broadcast spawning species (Ayre & Dufty, 1994; Hellberg, 1996; Maier et al., 2005; Nishikawa, Katoh & Sakai, 2003; Underwood et al., 2007). A. elisabethae is not an internal brooder, but developing embryos/larvae remain on the maternal colony surface for several days (Gutiérrez-Rodríguez & Lasker, 2004c) and are washed off of the colony as fully developed planulae. That should promote larval retention over long distance dispersal. Furthermore, the planulae are negatively buoyant, and 12% of planulae that were shaken from colonies landed on the substratum within a distance of 5 m of the natal colony in one set of observation in 1998 and 14% landed within 15 m of the natal colony during similar observations in 1999 (Gutiérrez-Rodríguez & Lasker, 2004c). Assuming a normal distribution, a value of σ of 83 m predicts that 5% of planulae would settle within 5 m and 14% settle within 15 m of the natal colony. Although unambiguous assignments were rare, recruits most commonly were assigned to either the local site or the adjacent site, again suggesting limited dispersal, albeit somewhat greater than that suggested by the IBD analysis.

Finding IBD among A. elisabethae populations was also dependent on sampling on an appropriate spatial scale. In this study, the proportion of variance explained by distance declined with the inclusion of more distant sites, which encompassed more complex seascapes. An earlier analysis (Gutiérrez-Rodríguez, Hannes & Lasker, 2005; Gutiérrez-Rodríguez & Lasker, 2004b), with a greater number of sites off of the LBB and using fewer markers identified a significant IBD pattern in the entire data set. However, the overall fit to the data was lower and the pattern among the LBB samples was not recognized. Subsequent reanalysis of those data indicates it was present in those data as well. Rosenberg et al. (2005) suggest that sampling intensity does not have a large effect on discerning structure. However, Selkoe & Toonen (2011) found sampling effects in analyses of IBD in a large survey of the literature on marine species. Similarly a number of authors have noted the importance of scale in differentiating IBD from other spatial patterns and report cases very similar to the A. elisabethae data (Meirmans, 2012; Nanninga et al., 2013; Polato et al., 2010; Selkoe et al., 2014; Skoglund et al., 2014; Treml et al., 2012; Tucker et al., 2014; Yasuda et al., 2014).

Seascape scale processes

Among many marine species simple geographic distance explains only a small portion of the variance in genetic distances, and many authors have adopted approaches emphasizing “seascapes” or “oceanographic distance” (Alberto et al., 2011; Galindo, Olson & Palumbi, 2006; Selkoe et al., 2010; White et al., 2010). Galindo, Olson & Palumbi (2006) found a good fit between the predictions of a biophysical model of larval dispersal with empirical data on the Caribbean coral Acropora cerviconis.White et al. (2010) found that consideration of ocean currents and multiple reproductive events could account for over 40% of the genetic differences of a whelk between sites in southern California, and Alberto et al. (2011) found that a measure of current driven connectivity that they called derived oceanographic distance was an excellent predictor of genetic differentiation. Similarly, both scleractinian corals (Baums, Paris & Cherubin, 2006; Foster et al., 2012) and reef fish (Taylor & Hellberg, 2006) exhibit structure predicted from oceanographic models. However, oceanographic distance alone does not offer a complete picture of gene flow. For instance, Foster et al. (2012) found that empirical and modeled estimates of gene flow do not correspond at smaller scales (e.g., Mesoamerican Barrier Reef System), and on a larger geographic and temporal scale, Crandall et al. (2014) concluded that contemporary dispersal predicted from currents was not adequate to explain gene flow of the sea star Linckia laevigata.

Among A. elisabethae populations in The Bahamas, the presence of IBD regressions with different slopes indicates that once sites away from the LBB were included, connectivity was not a simple function of distance. The degree of isolation between sites, while correlated with distance, probably integrates the effects of geography and current patterns in limiting dispersal. The pattern of IBD observed on the LBB extends across a 70 km open, deep water channel to EGG (Fig. 2B) but breaks down when more distant sites are included (Fig. 2A). Both observations illustrate the importance and potential limitations of seascape scale processes. Connectivity between EGG and the LBB may be driven by currents or eddies carrying larvae across the Northeast Providence Channel (the channel separating the eastern edges of the Little Bahama and Great Bahama Banks, Fig. 1). Alternatively, the population at EGG may have been founded recently and the similarity may reflect that historic event, not ongoing connectivity. In contrast, the other A. elisabethae sites appear to be isolated from the LBB, and isolated from each other. Although there was significant correlation between geographic and genetic distance when all of the sites are considered, closer inspection shows the data for the more distant sites (Exuma Sound (HC) and San Salvador (PR) in Fig. 2A) appear as groups of genetically distant sites with no pattern of IBD. Aside from PR (San Salvador) the other sites lie along the perimeter of the GBB, but the large expanse of shallow soft sediment environment of the GBB does not provide suitable habitat for stepping stone populations that could link the sites. Beyond the IBD pattern on the LBB and the distinctness of PR on San Salvador, the structure of the other sites is difficult to decipher. Models of current flow at the time of spawning in The Bahamas may provide further explanatory power in assessing the isolating mechanisms among populations.

Population structure of A. elisabethae in The Bahamas, is best described as sets of regional clusters. IBD can be found within areas such as the LBB when the distribution of reefs within those clusters and the sampling intensity is great enough. This illustrates how patterns of differentiation among populations are generated by a mix of IBD along reef tracts and isolating barriers occurring over larger distances and across oceanographic barriers. Nanninga et al. (2014) and Yasuda et al. (2014) reported similar cases in which IBD analysis of geographic subsets of sites generated slopes that differed from each other or from the slope of the combined data set. Polato et al. (2010) also report a similar pattern for the coral Porites lobata where a significant IBD pattern was identified as well as a Structure analysis that clearly divided the sites into two populations.

Sampling across multiple scales Ledoux et al. (2010) reported patterns of between region structure and within region IBD that are similar to in A. elisabethae. Similarly, the gorgonian Paramuricea clavata, which surface broods, exhibited distinct structure between regions across the Mediterranean and, as in our results, an IBD signal within some (Mokhtar-Jamai et al., 2013; Mokhtar-Jamai et al., 2011). The remarkably good fit of the A. elisabethae data to the IBD model compared to these other species is likely attributable to the linear pattern of suitable reef habitat along the Little Bahama Bank.

Selection and temporal variation

Genetic analyses of recruits estimates current patterns of gene flow and thus connectivity on an ecological/demographic time scale. The population structure of the A. elisabethae recruits closely matched that of the adults, and there was little indication of interanual differences in the pattern. In addition, a finer scale analysis of the same recruits Smilansky & Lasker (2014) found within site spatial autocorrelaton of the recruits in one of 3 years, which suggests that at times very local dispersal may dominate recruitment. These results would be expected in a species with limited dispersal ability such as A. elisabethae. Analysis of recruits over 2 to 6 years (CHSL and CHR, respectively) does not preclude the possibilty of unusual events but the similarity of those recruits to adults which live for decades (Goffredo & Lasker, 2005) indicates that the patterns of recruitment that we observed can account for the structure in the adult population. Similarly, the limited variation between years suggests that temporal variability in the delivery of recruits (Fodrie et al., 2011; Hogan et al., 2012) and “sweepstakes” effects (Carson et al., 2010; Pusack et al., 2014) were not occurring among these populations.

Analysis of recruits over multiple years also provides the opportunity to look for evidence of selection modulating the variability in settling recruits. Most of the recruits were collected 1 to 2 months after spawning. Selection against migrant genotypes could have occurred during that time period, but the lack of genetic difference between recruits and adult colonies suggests that maintenance of this pattern was more likely a function of local recruitment than selection limiting connectivity.

Management implications

The available data suggest that A. elisabethae larvae disperse over distances most commonly on the scale of 10² m. Some larvae may settle almost immediately (Gutiérrez-Rodríguez & Lasker, 2004c). The presence of A. elisabethae in Bermuda (Look et al., 1986a), demonstrates the distribution of dispersal distances has an extremely long tail, but both genetic (this study) and recruitment rate data suggest the majority of recruits are local in origin. Pinsky et al. (2012) suggest that populations will have closed demographics when the habitat is patchy on a scale >2σ, where σ is the standard deviation of the dispersal kernel. That describes the environment along the edges of both the Little Bahama and Great Bahama Bank. It suggests that population demographics is highly dependent on self-seeding. A. elisabethae is harvested in The Bahamas, with harvesters cropping colonies which are then allowed to recover. The movement of A. elisabethae larvae over distances of 10s and 100s of meters but not often further suggests that harvests or naturally occurring disturbances at those or greater scales will depress recruitment and slow recovery of populations. The slow recovery suggested by these patterns may also explain the patchy distribution of A. elisabethae where some sites have dense populations while other seemingly similar sites have few or none of the species. These differing densities may be a reflection of the history of the populations and not necessarily the suitability of the sites. Both the frequency and intensity of disturbances to marine systems and particularly coral reefs appears to be increasing (Gardner et al., 2005; Gardner et al., 2003). Understanding patterns of dispersal among marine species cannot reverse those trends but will be crucial to the development of management plans that can enhance recovery between disturbances. In the case of A. elisabethae low dispersal distances indicates that self-recruitment is important to the maintenance of populations. Incorporating the correct scaling in management plans will be critical to the maintenance of A. elisabethae and similar species.