Population genetic structure and connectivity of the seagrass Thalassia hemprichii in the Western Indian Ocean is influenced by predominant ocean currents

Abstract This study is the first large‐scale genetic population study of a widespread climax species of seagrass, Thalassia hemprichii, in the Western Indian Ocean (WIO). The aim was to understand genetic population structure and connectivity of T. hemprichii in relation to hydrodynamic features. We genotyped 205 individual seagrass shoots from 11 sites across the WIO, spanning over a distance of ~2,700 km, with twelve microsatellite markers. Seagrass shoots were sampled in Kenya, Tanzania (mainland and Zanzibar), Mozambique, and Madagascar: 4–26°S and 33–48°E. We assessed clonality and visualized genetic diversity and genetic population differentiation. We used Bayesian clustering approaches (TESS) to trace spatial ancestry of populations and used directional migration rates (DivMigrate) to identify sources of gene flow. We identified four genetically differentiated groups: (a) samples from the Zanzibar channel; (b) Mozambique; (c) Madagascar; and (d) the east coast of Zanzibar and Kenya. Significant pairwise population genetic differentiation was found among many sites. Isolation by distance was detected for the estimated magnitude of divergence (D EST), but the three predominant ocean current systems (i.e., East African Coastal Current, North East Madagascar Current, and the South Equatorial Current) also determine genetic connectivity and genetic structure. Directional migration rates indicate that Madagascar acts as an important source population. Overall, clonality was moderate to high with large differences among sampling sites, indicating relatively low, but spatially variable sexual reproduction rates. The strongest genetic break was identified for three sites in the Zanzibar channel. Although isolation by distance is present, this study suggests that the three regionally predominant ocean current systems (i.e., East African Coastal Current, North East Madagascar Current, and the South Equatorial Current) rather than distance determine genetic connectivity and structure of T. hemprichii in the WIO. If the goal is to maintain genetic connectivity of T. hemprichii within the WIO, conservation planning and implementation of marine protection should be considered at the regional scale—across national borders.


| INTRODUC TI ON
Understanding population genetic structure, long-distance dispersal and connectivity patterns of organisms facilitates conservation actions in spatially widespread coastal seascapes (Jones et al., 2009).
Population genetics may act as a powerful tool for resource management planners to understand the genetic connectivity between populations, which in turn may have implications for decisions regarding number, sizes, and locations of protected areas (Palumbi, 2003;Waycott et al., 2009). Furthermore, genetic studies can also provide estimates of genetic diversity, an important factor for an organism's adaptive capacity to survive in a changing environment (Smith & Bernatchez, 2008;Vandergast, Bohonak, Hathaway, Boys, & Fisher, 2008).
During the last decades, seagrass meadows have experienced a substantial decline (Orth et al., 2006;Waycott et al., 2009) and information regarding genetic structure and variability of meadows could provide an important link for coastal management actions. There are a few previous large-scale assessments of seagrass population genetics worldwide (see Arriesgado et al., 2015, Hernawan et al., 2016Jahnke et al., 2018, Sinclair et al., 2016Triest et al., 2018;van Dijk, Mellors, & Waycott, 2014, Wainwright, Arlyza, & Karl, 2018. While in the Western Indian Ocean (WIO) we are only aware of one local assessment of the seagrass Thalassodendron ciliatum from 2001 in southern Mozambique using RAPD markers (Bandeira & Nilsson, 2001) and one recent study about the seagrass Zostera capensis sampled along the South African coast, and in one bay in Mozambique and one bay in Kenya indicating the presence of two population clusters broadly corresponding to populations on the west and east coasts of Africa (Phair, Toonen, Knapp, & Heyden, 2019).
The seagrass Thalassia hemprichii (Ehrenberg) Ascherson is widely distributed in the Indo-Pacific (Green & Short, 2003). It is one of the most common seagrass species in the WIO (Bandeira & Björk, 2001;Gullström et al., 2002), a biogeographic sub-region of the Indian Ocean stretching on a latitudinal scale from Somalia to the east coast of South Africa (Obura, 2012). T. hemprichii reproduces sexually by seeds and asexually by rhizome growth. The extent of each component of reproduction has important effects on local population demographics, dispersal, biogeography, and genetic diversity. Seed banks, in terms of long-term survival of buried seeds, play an important role in the persistence of some seagrass species, but seem to be absent in T. hemprichii (Rollon, Vermaat, & Nacorda, 2003). In terms of long-distance dispersal, positively buoyant shoots with attached rhizomes or seedling have the highest potential for long-distance dispersal in Thalassia spp. (Kendrick et al., 2012;Wu, Chen, & Soong, 2016). In the west Pacific, adult plants of T. hemprichii have been shown to be able to float for months and still remain alive and potentially able to colonize new areas (Wu et al., 2016). Seeds seem to sink within 24 hr and fruits, while having a potential to float for about a month and at distances of dozens to hundreds of kilometers, do not seem to contain seeds that can germinate after long-time floatation (van Dijk, Tussenbroek, Jiménez-Durán, Márquez-Guzmán, & Ouborg, 2009;Lacap, Vermaat, Rollon, & Nacorda, 2002;Wu et al., 2016). The fruits and seeds do not seem to survive the passage through the birds' digestive tract (Wu et al., 2016). Therefore, long-distance dispersal seems to be primarily dependent on buoyant shoots (and to a lesser degree on fruits) that move passively by currents (Thiel & Gutow, 2005). Hence, oceanographic features are expected to strongly shape the directionality of dispersal and/or to act as barriers for T. hemprichii.
Major currents in the WIO region are driven by the east-west flow of the South Equatorial Current (SEC) that carries water toward Madagascar, where it splits up into the North East Madagascar Current (NEMC) continuing toward the southern coast of Tanzania and the northern coast of Mozambique (Obura et al., 2018;Richmond, 2002) and the South East Madagascar Current (SEMC) continuing south along the east coast of Madagascar. Where the NEMC reaches the coast of southern Tanzania and northern Mozambique, this major current system splits into two major currents, that is, the East African Coastal Current (EACC) that flows northwards along the coasts of Tanzania and Kenya and the Mozambique current system that flows south along the coast of Mozambique creating complex eddies in the Mozambique Channel (Figure 1; Benny, 2002;Obura et al., 2018;Richmond, 2002).

K E Y W O R D S
coastal conservation, connectivity, dispersal, gene flow, genetic structure, microsatellite, ocean current, population genetics, seagrass, Western Indian Ocean genetic sub-region of the Indian Ocean (Ridgway & Sampayo, 2005). This makes regional population genetic studies highly valuable for management initiatives in the WIO (Ridgway & Sampayo, 2005).
The main aim of the current study was to reveal the genetic population structure, dispersal patterns, and connectivity of the common and important habitat forming seagrass T. hemprichii in the WIO region, with focus on the coasts of Kenya, Tanzania (mainland and Zanzibar), Mozambique, and Madagascar. Since dispersal of both shoots and fruits of T. hemprichii mainly occurs by passive rafting, we hypothesized that genetic population structure in the region is shaped by the South Equatorial, the North East Madagascar, and the East African Coastal currents. The findings are used to discuss how to improve seagrass management in the WIO region to sustain healthy seagrass populations with high resilience to environmental change.

| Study area and field sampling
The WIO covers about 40 degrees of latitude and hosts a high level of marine biodiversity (from Somalia in the north to eastern South Africa in the south). The climate and pattern (and strength) of currents in the WIO are complex and strongly influenced by the monsoonal circulation ( Figure 1). Throughout the SE monsoon season (March to October), the EACC is speeded up by the southeasterly trade winds and consequently, the current speed reaches ca. 1.5-2 ms −1 . During the NE monsoon season (October-March), the current is slowed down by the north-easterly trade winds and as a consequence the current speed is lowered to about 0.5 ms −1 (McClanahan, 1988;Mahongo & Shaghude, 2014).   Figure 1).
One leaf per shoot and one piece of rhizome were cleaned of epiphytes and dried, and stored in silica crystals or frozen in 70% ethanol until DNA extraction.

| Scoring and data quality checks
Fragment size analysis was performed with undiluted PCR products using an Applied Biosystems 3730 DNA Analyser with a 350 ROX internal size standard added to each well. The fragments were scored with PeakScanner ® . Samples with an unclear allele pattern were reamplified and re-genotyped. We did not succeed in amplifying all individuals at all loci, and individuals with missing data were removed.
When working with partially clonal organisms, it is important to distinguish between unique genotypes and clonal replicates, and we used RClone (Bailleul, Stoeckel, & Arnaud-Haond, 2016) in R 3.3.1 (R Development Core Team, 2014) to identify multilocus genotypes (MLGs), that is clones. As most downstream analyses assume microsatellite loci to be independent from each other (i.e., not in linkage) and to be in Hardy-Weinberg equilibrium (HWE), we tested for HWE at each locus and across all loci in each population with Genepop 4.2 (Raymond & Rousset, 1995) using 100 batches and 1,000 iterations per batch, and applying Bonferroni corrections, and for Linkage Disequilibrium (LD), a modification of the standard methods, described by Agapow and Burt (2001) in the R package poppr (Kamvar, Tabima, & Grünwald, 2014), was applied using 1,000 permutations.
Loci with null alleles (i.e., not amplified allele copies) need to be removed from the dataset since they cause genotyping errors when one of the allele copies fails to be amplified by the PCR, which leads to missing genotypes and the incorrect assignment of homozygotes.
To detect null alleles, we used the software MicroDrop (Wang & Rosenberg, 2012) with 10,000 permutations and 100 replicates.
Furthermore, molecular markers used to assess connectivity should not be under selection, and therefore it is recommended to always test for outlier loci before estimating population genetic parameters To test if any of our loci deviate significantly from expectation under neutrality, we ran BayeScan (Foll & Gaggiotti, 2008) with default settings and used the R script provided by Foll and Gaggiotti (2008) to identify any loci showing signs of selection.

| Genetic diversity and population differentiation
To estimate genetic diversity (genetic variation within each sampling site) and population differentiation (genetic difference between sampling sites) several calculations were made. Genotypic richness (R), a measure of the amount of clonal growth, was calculated with the formula (MLG − 1)/(N − 1), which relates the number of MLGs to the number of ramets (N) (Dorken & Eckert, 2001 (Weir & Cockerham, 1984) and the unbiased estimator of Jost D EST (Jost, 2008) and G ST ′ (Hedrick, 2005)  In order to detect spatial genetic structure in the WIO, we used Bayesian clustering algorithms implemented in STRUCTURE (Pritchard, Stephens, & Donnelly, 2000) and TESS 2.3 (Chen, Durand, Forbes, & François, 2007). Both programs perform Markov Chain Monte Carlo simulations with the aim to find the most probable number of genetic clusters under Hardy-Weinberg equilibrium.
As we had geographic information only at the population level, we used TESS to calculate slightly adapted geographic coordinates for each individual. We then run TESS with the CAR admixture model, which assumes spatial autocorrelation of genetic differentiation, using the default value of 0.6 for the strength of the autocorrelation.
We run TESS for an assumed number of clusters, K max = 2-12 using a burn-in of 10,000 sweeps followed by 25,000 sweeps, with 100 independent runs conducted for each K max . We then used the average deviance information criterion for each value of K max to evaluate the most likely number of genetic clusters. For visualizing clusters and for post-processing of TESS outputs, we used pophelper (Francis, 2016) with CLUMPP 1.1.2 (Jakobsson & Rosenberg, 2007) in R 3.2.2.

| Directional migration rates and bottlenecks
To calculate directional migration rates among sampling sites, we used divMigrate-online (https ://popgen.shiny apps.io/divMi grateonlin e/) based on the genetic distance measures of G ST (Sundqvist, Keenan, Zackrisson, Prodohl, & Kleinhans, 2016). To identify firstgeneration migrants (individuals with different allele frequencies than the population where they were sampled), we also used an assignment test implemented in GENECLASS2 (Piry et al., 2004) with the exclusion method. We investigated local population dynamics and the potential presence of genetic bottlenecks with the program Bottleneck 1.2.02 (Cornuet & Luikart, 1996) assuming a step-wise mutation model (SMM) and the two-phased model of mutation (TPM) and using 10,000 replications. We only conducted this test at sites with more than 10 MLGs.

| Isolation by distance
Isolation by distance (IBD) is the relationship of genetic similarity between sites over a geographic distance. The model tests the assumption that the further apart two sites are from each other geographically, the more genetically different they should be. To test for IBD, we correlated the genetic distance matrices with "sea distance" (geographic distance among sampling sites without crossing land).
This measure was calculated with the R package marmap (Pante & Simon-Bouhet, 2013). Subsequently, we carried out a Mantel test using the R package ncf (Bjornstad, 2009) in R 3.2.2 and resampling of the matrices was performed (100,000 times).

| Genetic data quality control
Clonality varied among sites and we identified 6-19 genotypes per population (

| Genetic diversity and differentiation
Genetic diversity was generally low.  (Weir & Cockerham, 1984) ranged between −0.05 and 0.64, and two-thirds of the pair-
The additional assignment test identified seven first-or secondgeneration migrants, of which two could be assigned to other sampling sites ( Table 2). The immigrant at the site ZP in Mozambique is predicted to come from either TZF or TZM (Zanzibar), as also suggested by the divMigrate analysis ( Figure 4). The immigrant at F I G U R E 2 The PCA analysis of the 11 Thalassia hemprichii populations in the Western Indian Ocean shows a horseshoe pattern of differentiation on the first two axes. Samples coming from the four sites in Zanzibar (starting with TZ) are colored in shades of red, the Tanzanian mainland site (TM) in black, the three sites from Mozambique (starting with Z) in shades of blue, the two sites from Madagascar (starting with M) in shades of green, and the Kenyan sampling site (KM) in purple. Site acronyms as in Table 1 and Figure 1; the inset plot shows the eigenvalues of the first two principal components, which were used in this analysis F I G U R E 3 TESS clustering analysis of the 11 Thalassia hemprichii populations in the Western Indian Ocean. The TESS plot is shown for the most likely number of clusters (according to the deviance information criterion of averaged runs) and plots for other K max can be found in Appendix S6. Within each plot, each vertical bar represents an individual belonging to the sampling location indicated under the x-axis, clusters are color coded, and the y-axis of each plot shows the proportion of the genotype belonging to each cluster. Site acronyms as in Table 1

The bottleneck tests indicated that TZC on Zanzibar and ZIS in
Mozambique might have experienced a recent population decline (pvalue < 0.05 for TPM and SMM) and TZF on Zanzibar exhibited signs of a recent expansion (p-value < 0.05 for TPM and SMM).

| Isolation by distance
"Sea distance" ranged from ~15 to ~2,700 km. The two Madagascan sites are most geographically close, while the samples from Inhaca, Mozambique (ZIN and ZIS), and the Kenyan site (KM) are furthest apart. There was no evidence for a significant correlation between sea distance and the genetic distance measure F ST (p = 0.17, r = 0.14) nor G' ST (p = 0.14, r = 0.17), but sea distance is significantly correlated with Jost's D EST (p = 0.03, r = 0.37). Although most pairwise F ST values are significant, D EST is arguably the most meaningful measure for our dataset given that the observed heterozygosity is very low.

| D ISCUSS I ON
The overall aim of this study was to understand the large-scale (1000s km) population genetic pattern of the seagrass T. hemprichii in the WIO region for regional conservation purposes. The study contributes to the knowledge of population genetic patterns of important species in the region, but also where in fact little is known on population genetic patterns in general. Genotyping of the seagrass T. hemprichii was successful, and we found four distinguished genetic clusters separated geographically on various spatial scales (10s-1000s km). Interestingly, a clear genetic cluster co-occurs with a local system of eddies between the mainland of Tanzania and Zanzibar. Gene flow was directional and strongest from Madagascar toward the coasts of Kenya and western Zanzibar, spatially coinciding with the NEMC and the "continuing" EACC (Zavala-Garay et al., 2015). Our findings indicate that the genetic structure of T. hemprichii in the WIO is influenced by large oceanic currents (SEC, NEMC, and EACC), as well as by local hydrodynamic patterns (Zanzibar channel). The study, therefore, emphasizes the necessity to consider different spatial scales to understand metapopulation dynamics and identify source populations of a habitat providing plant species.

| Genetic structure and connectivity of T. hemprichii in the WIO
We show significant genetic structure among most of the sampled T. hemprichii meadows and find a relation between genetic differentiation (D EST ) and increasing geographic distance. Nevertheless, large-scale currents in the area seem to play an important role for the observed genetic differentiation and the four genetic clusters suggested by the TESS analysis conform well with the major oceanographic currents in the WIO. We hypothesized that the northbound TA B L E 2 Assignment test of Thalassia hemprichii at 11 locations in the West Indian Ocean Note: For each site (acronyms as in Table 1), the sampling sites are presented in rows and the percentage of individuals that is assigned to their own sampling site is shown (% assigned own site). A maximum of two migrants were found at one sampling site and the putative origin of each migrant is shown. When more than one location showed likelihoods above the threshold of 0.1, the most probable as well as other probable sources are shown. We calculated the probability that an individual belongs to the population from which it was sampled using a partially Bayesian criterion (Rannala & Mountain, 1997) and compared the likelihood of exclusion of an individual to a distribution of likelihoods of 10,000 simulated genotypes in order to define a statistical threshold (Paetkau, Slade, Burden, & Estoup, 2004;Underwood, Smith, Oppen, & Gilmour, 2007) with a type I error of 0.01. We excluded an individual from its sampling site when the probability for exclusion was above 95% and then excluded the migrants from the dataset, which served as the reference to which migrants were assigned. We assigned migrants to another sampled population when the probability was p ≥ 10% for this other population (Underwood et al., 2007).  (Souter et al., 2009;van der Ven et al., 2016), and therefore it seems likely that the sites in the Zanzibar channel are subject to local oceanographic conditions that result in propagule retention and oceanographic isolation.
A higher admixture was found for the sites in Mozambique, Zanzibar, and Kenya, which is in concordance with other genetic studies in the area on crabs (Fratini, Ragionieri, & Cannicci, 2010), corals (van der Ven et al., 2016, and fish (Visram et al., 2010).

| Contribution of clones and heterozygote deficits
Asexual reproduction is common in seagrasses and has knock-on effects on genetic variation of populations. We detected high genotypic richness at most sites, despite distances of 10-150 m between individuals, which suggests that both sexual and asex-  (Hernawan et al., 2016). The combination of low marker polymorphism and high clonality limits the power of the genetic inference analyses performed here (Hale, Burg, & Steeves, 2012;Ryman et al., 2006) and give evidence for the importance of developing markers that are more informative in the WIO. In particular, the combination of high clonality and low marker polymorphism found in this study results in a high likelihood that a large proportion of alleles are not detected (Hale et al., 2012), which has impacts on all F ST related measures. As we observed very high values from the F ST measures and low polymorphism, we used the approach suggested by Wang (2015) to test whether our estimates are substantially affected by mutations. Results (not shown) revealed no significant negative relationship between the two measures, indicating that results are reliable (Wang, 2015).
A previous study on the clown fish Amphiprion akallopisos showed higher genetic diversity in the East Indian Ocean (EIO) compared to the WIO, suggesting that this fish species originated in the EIO (Huyghe & Kochzius, 2017). Similarly, we hypothesize that the origin of T. hemprichii is in the EIO and Pacific and low genetic diversity may be explained by the fact that our sampling was performed at the edge of its distributional range. Furthermore, a number of studies have shown limited gene flow between the WIO and the EIO and Pacific, while others have shown low genetic diversity in the WIO (reviewed in Ridgway & Sampayo, 2005).

| Management implications and conclusion
This study is the first large-scale genetic study of a seagrass in the WIO and results indicate that dispersal of T. hemprichii is to a large extent influenced by predominant currents in the area. Nevertheless, traditional isolation by distance models-where genetic differentiation increases with geographical distance-also seems to play a role at the assessed spatial scales in the WIO, as do local hydrodynamic features. For future conservation planning, the information on genetic clusters and breaks, predominant directions of gene flow and genetic diversity should be used to ensure that connectivity is maintained in the WIO. A regional scale approach to seagrass conservation is recommended, including effective cross-boundary management.

ACK N OWLED G M ENTS
We thank Jeanine L. Olsen for the use of her former laboratory and  -1288SWE-2010-194;SWE-2012-086). In Tanzania and Madagascar, in-country authors and staff, associated with in-country institutions, handled the seagrass. Most of the seagrass collection was carried out before the Nagoya protocol entered into force. Nevertheless, all data are openly available for fair and equitable sharing of benefits. We thank the anonymous reviewers for constructive comments on an earlier version of the manuscript.

CO N FLI C T O F I NTE R E S T
None declared.