Genetic population structure and demography of an apex predator, the tiger shark Galeocerdo cuvier

Abstract Population genetics has been increasingly applied to study large sharks over the last decade. Whilst large shark species are often difficult to study with direct methods, improved knowledge is needed for both population management and conservation, especially for species vulnerable to anthropogenic and climatic impacts. The tiger shark, Galeocerdo cuvier, is an apex predator known to play important direct and indirect roles in tropical and subtropical marine ecosystems. While the global and Indo‐West Pacific population genetic structure of this species has recently been investigated, questions remain over population structure and demographic history within the western Indian (WIO) and within the western Pacific Oceans (WPO). To address the knowledge gap in tiger shark regional population structures, the genetic diversity of 286 individuals sampled in seven localities was investigated using 27 microsatellite loci and three mitochondrial genes (CR,COI, and cytb). A weak genetic differentiation was observed between the WIO and the WPO, suggesting high genetic connectivity. This result agrees with previous studies and highlights the importance of the pelagic behavior of this species to ensure gene flow. Using approximate Bayesian computation to couple information from both nuclear and mitochondrial markers, evidence of a recent bottleneck in the Holocene (2,000–3,000 years ago) was found, which is the most probable cause for the low genetic diversity observed. A contemporary effective population size as low as 111 [43,369] was estimated during the bottleneck. Together, these results indicate low genetic diversity that may reflect a vulnerable population sensitive to regional pressures. Conservation measures are thus needed to protect a species that is classified as Near Threatened.


| INTRODUC TI ON
Study of large sharks, including the tiger shark Galeocerdo cuvier, the great white shark Carcharodon carcharias, and the whale shark Rhincodon typus, is challenging as these species spend substantial periods of their lifetime in open ocean waters. Consequently, data concerning basic aspects of their biology, such as migration patterns and population structure, are limited (Conrath, Musick, Carrier, & Heithaus, 2012;Musick, 2010). Understanding the ecological role of large sharks in marine ecosystems also remains incomplete. However, as apex predators they are considered to exercise important functions in marine food webs via top-down processes (Dudley & Simpfendorfer, 2006;Ferretti, Worm, Britten, Heithaus, & Lotze, 2010;Myers, Baum, Shepherd, Powers, & Peterson, 2007). Large sharks typically present classically K-selected life histories, with slow growth rate, late maturity, and low fecundity (Musick, Burgess, Cailliet, Camhi, & Fordham, 2000). This renders them vulnerable to overexploitation with low rebound potentials limiting their recovery (Campana & Ferretti, 2016;Cortés, 2002;Dudley & Simpfendorfer, 2006;Ferretti et al., 2010;Myers & Worm, 2005;Worm et al., 2013). Certain species, including the white, tiger, and bull sharks, are also responsible for the majority of human-shark conflicts, complicating conservation and management actions. Consequently, continuing to build on our current understanding of the biology and the ecology of these large sharks is needed to facilitate appropriate management .
The tiger shark is a large (up to 5.5 m long) Carcharhinid (Meyer et al., 2014) with a circumglobal distribution in tropical and subtropical waters (Compagno, 1984(Compagno, , 1990. It is a potential keystone species in marine ecosystems through predation or by inducing behavioral modifications of its prey, and thus indirectly modifying primary producer community structure, biomass, and nutrient composition (Burkholder, Heithaus, Fourqurean, Wirsing, & Dill, 2013;Heithaus, Frid, Wirsing, & Worm, 2008;. The tiger shark is listed as globally "Near Threatened" by the International Union for Conservation of Nature (IUCN) Red List of Threatened Species and is primarily threatened by fisheries exploitation (Simpfendorfer, 2009;Temple et al., 2018). Clarke et al. (2006) estimated that ~400,000-500,000 tiger sharks are caught annually for the shark fin trade globally. Furthermore, this species is currently targeted by shark control programmes in the Indo-Pacific: South Africa (Cliff & Dudley, 1991;Dudley, 1997;Sumpton, Taylor, Gribble, McPherson, & Ham, 2011) and Australia (Holmes et al., 2012;Reid & Krogh, 1992;Simpfendorfer, 1992), and formerly in Hawaii (Wetherbee, Lowe, & Crow, 1994). It is also reported as bycatch in pelagic fisheries, in the Western Pacific (Polovina & Lau, 1993) and in the southern (Afonso & Hazin, 2014) and northwestern Atlantic (Baum et al., 2003). Trends in long-term catch and catch rates are difficult to obtain for sharks not specifically targeted by fisheries. Nevertheless, long-standing control programmes as well as logbooks from the longline fishing fleets do provide long-term data using standardized fishing methods, which may be used to assess catch per unit effort trends over time.
Tiger shark catch rates in control programmes appear to be increasing in KwaZulu-Natal, South Africa (Dudley & Simpfendorfer, 2006), while contrasting declines were observed in Queensland (Holmes et al., 2012), New South Wales, Australia (Reid, Robbins, & Peddemors, 2011), and in commercial fisheries in the northern Atlantic (Baum et al., 2003;Myers et al., 2007). These differences in catch rates indicate regional variation in population trends of tiger sharks, but overall evidence supports declining populations.
Control programme catch rates coupled with mark-recapture and telemetry also provide information on habitat use patterns.
Data from pelagic longline fisheries have highlighted the importance of the pelagic realm to a species originally described as coastal (Domingo et al., 2016;Polovina & Lau, 1993). More recently, tiger shark individuals have been recorded moving over several thousands of kilometers distances (Ferreira et al., 2015;Hammerschlag, Gallagher, Wester, Luo, & Ault, 2012;Holmes et al., 2014;Lea et al., 2015;Werry et al., 2014), including crossing ocean basins (Afonso, Garla, & Hazin, 2017;Kohler, Casey, & Turner, 1998;Kohler & Turner, 2001). Equally, tracking studies have also revealed strong residency patterns, with some individuals maintaining large but defined home ranges and returning to specific locations on a regular basis (Ferreira et al., 2015;Fitzpatrick et al., 2012;Heithaus, 2001;Holland, Wetherbee, Lowe, & Meyer, 1999;Lowe, Wetherbee, & Meyer, 2006). These patterns seem to be linked not only to intrinsic states such as foraging strategies (Heithaus, Hamilton, Wirsing, & Dill, 2006;Holland et al., 1999;Meyer, Clark, Papastamatiou, Whitney, & Holland, 2009;Meyer, Papastamatiou, & Holland, 2010;Papastamatiou et al., 2011) and sex Papastamatiou et al., 2013;Sulikowski et al., 2016) but also to extrinsic drivers, notably prey abundance (Lowe, Wetherbee, Crow, & Tester, 1996) and water temperature (Heithaus, 2001;Holmes et al., 2014;. While these movement studies highlight the complex migration patterns and habitat use of tiger sharks, the implications of these patterns on population connectivity and structure have only begun to be considered recently using molecular markers. The first study investigating the population genetic structure of the tiger shark was primarily designed for species delimitation examining 29 samples collected throughout its range using the mitochondrial NADH dehydrogenase subunit (ND2) gene (Naylor et al., 2012). Two monophyletic clades were identified, one in the Atlantic and the other in the Indo-Pacific, with no shared haplotype, suggesting the presence of two subspecies. This hypothesis was subsequently refuted by Bernard et al. (2016), using a larger sample of 380 individuals from several sampling sites across the three ocean basins, 10 microsatellite loci and two mitochondrial genes, the control region (CR) and the cytochrome oxidase c subunit I (COI). Bernard et al. (2016) highlighted long-term genetic isolation between tiger shark populations of the Atlantic and the Indo-Pacific, but with shared mitochondrial haplotypes, which is inconsistent with the two subspecies hypothesis. Furthermore, samples from Hawaii appeared genetically differentiated from all other locations, perhaps due to more restrictive movement behaviors and greater residency exhibited by sharks from this area (Meyer et al., 2010;Papastamatiou et al., 2013). Using mitochondrial data, Bernard et al. (2016) also identified a larger degree of genetic differentiation in the Indo-Pacific than observed with microsatellite data. Recorded differences between the west and east coast of Australia were hypothesized to be a result of a stronger matrilineal structure due to sex-biased dispersal and female philopatry. An additional study focused on the population structure of tiger sharks across the eastern Indian Ocean and the Pacific, with 355 samples collected primarily around Australia and Hawaii, with eight samples from Brazil in the southwest Atlantic as an outgroup (Holmes et al., 2017). Using nine microsatellite loci, Holmes et al. (2017) confirmed genetic isolation between the western Atlantic and the Indo-Pacific.
Further, the high connectivity around Australia found by Bernard et al. (2016) was reflected in satellite tracking data (Afonso et al., 2017;Heithaus et al., 2007;Kohler et al., 1998). Nevertheless, contrary to the findings of Bernard et al. (2016), Holmes et al. (2017) found no genetic differentiation between Hawaiian and Australian populations.
In the present study, we used samples from the western Indian Ocean (233 individuals from four locations), the eastern Indian Ocean (nine individuals from one location), and the western Pacific (33 individuals from two locations), and 27 microsatellite loci and three mitochondrial genes (CR, COI, and cytb) to further investigate population structure and demographic parameters in the tiger shark.
As we used samples in common with Holmes et al. (2017) and the nine microsatellite loci they used, we were able to combine samples to obtain a more precise picture of the population genetic dynamics of the tiger shark in the Indo-Pacific. In addition, CR sequences obtained by Bernard et al. (2016) were used in conjunction with our data, to further investigate their proposed mitochondrial structure in these oceans.

| Sampling
Samples were collected in four locations in the western Indian Ocean (Zanzibar, ZAN: n = 8; South Africa, SAF: n = 34; the Seychelles, SEY: n = 24; Reunion Island, RUN: n = 167), and in the eastern Indian Ocean, along the Australian west coast (AUS1: n = 9).
Samples were collected in the western Pacific from the northeast coast of Australia (Queensland, AUS2: n = 10) and in New Caledonia (NCA: n = 23; Figure 1). Samples came from individuals caught by fishermen or shark control programmes (fin clips or muscle tissue) and from scientific projects (i.e., biopsies) and were preserved in 90% ethanol.

| Laboratory procedures
Genomic DNA was extracted using Qiagen DNeasy Blood & Tissue kit (Qiagen, Hilden, Germany) following manufacturer instructions.
Of the nine loci developed by Bernard et al. (2015), we did not keep TGR233 as we found it difficult to read. Of the loci developed by Mendes et al. (2016) Holmes et al. (2017) with the nine microsatellite loci developed by Bernard et al. (2015). Blue indicates samples sequenced at the control region by Bernard et al. (2016) to read with several failed amplifications, or TIG25 as it was monomorphic throughout our samples. We thus kept 27 microsatellite loci for further analyses.
Microsatellite genotypes generated by Holmes et al. (2017) and Bernard et al. (2016) were also added to the present study ( Figure 1; Table 1). Individuals from Australia genotyped in this study are the same as those in Holmes et al. (2017), which were genotyped with the nine microsatellite loci developed by Bernard et al. (2015). For these individuals and the eight microsatellites in common between the present study and Holmes et al. (2017), the genotypes were compared and allele lengths calibrated. To ensure the allele frequency bins were uniform between the studies at each locus, alleles frequencies were plotted and compared for each sampling location (Appendix S2). We thus added genotypes of all the individuals from Holmes et al. (2017) (Figure 1, Table 1) enlarging the geographic coverage of our sampling: northern and TA B L E 1 Summary of Galeocerdo cuvier sampling locations and number of individuals from this study, Holmes et al. (2017) and Bernard et al. (2016), as well as the molecular markers used in each study and the different datasets analyzed in the present study We also expanded the number of individuals for some locations in common: western (AUS1) and northeastern (AUS2) Australian coasts, and New Caledonia ( Figure 1, Table 1).
Moreover, for DAPC analyses (see below) we also used the microsatellite genotypes generated by Bernard et al. (2016), retaining only the eight microsatellites in common with our study and Holmes et al. (2017). However, we could not analyze them together with our genotypes as we did not have in our possession samples from individuals in common to calibrate the electrophoregrams.

| Genetic diversity analyses
The mean allelic richness A r and the mean private allelic richness A rp were calculated using a rarefaction method, as implemented in HP-rare v.1.0 (Kalinowski, 2005). This method accounts for differences in sample size by standardizing A r and A rp values across sampled locations by resampling the lowest number of genotypes available (i.e., 16 haploid gene copies or eight diploid genotypes in Zanzibar) in each location.
Mitochondrial sequences were quality checked and aligned performed assuming a HKY85 model of substitution as the latter was shown to best fit the data (see Section 3). Rate of variation among sites was modeled with a discrete gamma distribution with four rate categories. We assumed an uncorrelated lognormal relaxed clock to account for rate variation among lineages. To minimize prior assumptions about demographic history, we adopted an extended Bayesian skyline plot (EBSP) approach in order to integrate data over different coalescent histories. Evolutionary model parameters were then estimated, with samples drawn from the posterior every 10 5 MCMC steps over a total of 10 8 steps from five independent runs.

| Population genetic structure
Two complementary clustering methods were used to investigate population structure in the tiger shark. First, Bayesian clustering analyses were performed using Structure v.2.3.4 (Falush, Stephens, & Pritchard, 2003;Pritchard, Stephens, & Donnelly, 2000). For any given number of clusters (K) between 1 and 10, individual assignment probabilities to each cluster were determined so as to minimize departures from HWE within clusters and maximize LD among clusters. Two analyses were performed, with and without the LOCPRIOR model, which uses prior sampling location information in the Bayesian clustering, to allow detection of weaker genetic population structure (Hubisz, Falush, Stephens, & Pritchard, 2009). Conditions were set to 10 6 chain length after a burn-in of 5 × 10 5 and 10 chains were run for each K, assuming correlated allele frequencies and the admixture model. The analysis was also performed using the alternative ancestry prior as suggested in Wang (2017) Assessing population differentiation between pairs of sampling locations, F ST (Weir & Cockerham, 1984) and D est (Jost, 2008) were estimated for the microsatellite data with arlequin  (Slatkin, 1995) was estimated using arlequin v.3.5.1.2 (Excoffier & Lischer, 2010). Significance of pairwise population differentiation indices was tested using 10,000 permutations.

| Population demography
To test for departures from a constant population size (Ramos-Onsins & Rozas, 2000), the summary statistics Tajima's D (Tajima, 1989) and Fu's F S (Fu, 1997) were estimated from the CR-COI-cytb dataset with arlequin v.3.5.1.2 (Excoffier & Lischer, 2010), with significance tested implementing 10 5 simulated samples. The time was measured backward in generations before present. N 0 , the actual effective population size; N 1 , the ancestral effective population size; N b , the effective population size during a bottleneck; N e , the effective population size during an expansion; t 1 , beginning of decrease or expansion for Scenarios 1 and 3; t 2 , beginning of decrease or expansion for Scenarios 2 and 4; t, beginning of the expansion or bottleneck period for Scenarios 5 and 6 (less than 500 generations) expansion (N 0 > N 1 ), (Scenario 4) a more ancient (between 10 3 and 5 × 10 5 generations in the past) expansion, (Scenario 5) an expansion followed by a decrease (N e &gt; N 0 ,N 1 ), (Scenario 6) a bottleneck (N b <N 0 ,N 1 ), and (Scenario 7) a constant effective population size (N 0 = N 1 ). The parameters t 1 (t 1 < 500 generations i.e., ≈5,000 years) and t 2 (t 2 > 500 generations) were chosen to reflect relatively recent events that may be linked to anthropogenic factors, or more ancient events, such as glacial/interglacial transitions. For Scenarios 5 and 6, the end of the expansion/bottleneck was set at five generations in the past, which approximately corresponds to the ban on commercial exploitation of the tiger shark in Reunion Island (in 1999). Generation time was supposed to be around 7-10 years (Branstetter, Musick, & Colvocoresses, 1987;Holmes et al., 2015;Kneebone, Natanson, Andrews, & Howell, 2008;Wintner & Dudley, 2000).
To run this analysis, we considered all the individuals from Reunion Island, as it was the location with the highest number of individuals and may be the most representative of the genetic diversity of the whole population. It has indeed been shown that pooling individuals from different sampling locations, even with nonsignificant pairwise differentiation values may bias results (Lombaert et al., 2014). Nevertheless, the analysis was also run using only samples from New Caledonia (n = 23), to consolidate the results.
For each scenario, 10 6 simulated datasets were run. To select the best scenario, posterior probabilities were computed via logistic regression on the 1% of simulated datasets closest to the empirical data (Cornuet et al., 2008). Summary statistics were transformed by linear discrimination analysis prior to logistic regression to reduce correlation among explanatory variables and provide conservative estimates of scenario discrimination (Estoup et al., 2012). Posterior distributions of all parameters were then estimated from the selected model, based on the 1% of simulated datasets closest to the empirical data. More details on the ABC analysis are provided in Appendix S3.

| 27-msat dataset
Null alleles were detected for several loci in several sampling loca-

| Mitochondrial diversity indices for the three mitochondrial genes
We obtained sequences of 869 bp for CR, 652 pb for COI and 931 bp for cytb.
Using only our samples, the three genes resolved 9, 4, 14 haplotypes with 7, 3, and 12 polymorphic sites, respectively. Total haplotype diversity (h) was lower for COI (0.09 ± 00) than for cytb (0.47 ± 0.00) and CR (0.48 ± 0.00). Similar results were observed for the total nucleotide diversity (π), varying from 0.00014 ± 0.00000 for COI to 0.00063 ± 0.00001 and 0.00068 ± 0.00000 for cytb and CR, respectively. Overall, variations of haplotype and nucleotide diversities were not constant across locations and across genes, with lowest values estimated for Reunion Island (h = 0.34 ± 0.01 and π = 0.00047 ± 0.00005) and Zanzibar (h = 0.00 ± 0.00 and π = 0.00000 ± 0.00000) and highest values for AUS1 (Western Australia; h = 0.78 ± 0.44 and π = 0.00115 ± 0.00032) for CR and TA B L E 3 Summary statistics for each sampling location of Galeocerdo cuvier for (a) the 2,452 bp CR-COI-cytb dataset and (b) the CR dataset: N S , number of individuals sequenced in this study; N Bernard , number of individuals sequenced in Bernard et al. (2016)

| Mitochondrial diversity on the CR dataset
Within the CR dataset, we analyzed 538 CR sequences and resolved 25 haplotypes and 16 polymorphic sites, with an overall haplotype diversity h of 0.74 ± 0.00 and a nucleotide diversity π of 0.00280 ± 0.00000 (Table 3b)

| Population genetic structure
Performing Bayesian clustering analyses without the LOCPRIOR model, no distinct genetic clusters were retrieved in either the 27msat or the 8-msat datasets (Figure 4a, c).    Table 5).  Concerning the ABC analyses (Appendix S10), observed summary statistics for all scenarios fell within the distribution of simulated summary statistics, suggesting adequate choice of prior distributions (Appendix S10a). Scenario 6 (population decrease at time t from N 1 to N b , followed by an expansion from N b to N 0 five generations ago; Figure 2) presented the highest posterior probability based both on the logistic regression-based estimates and the direct estimate of posterior probability (Appendices S10b and S10c).

| Population demography
Other scenarios received no statistical support. Furthermore, posterior error rates were relatively low, with values of 0.342 and 0.294 using the direct and logistic approaches, respectively. We thus estimated parameter values using data simulated under Scenario 6 ( to approximately 2,000 -3,000 years (RMSE = 1.347). However, the actual effective population size N 0 was less precisely inferred (Table 6; Appendix S11).
Similar results were found when performing the analysis including only individuals from New Caledonia (Appendices S12, S13, and S14).
Nevertheless, achieving an extensive sampling is difficult, especially for elusive species, while elaborating numerous molecular markers remains expensive. Thus, researchers continually strive to find a balance between both strategies but ultimately have to settle for what they can afford. Hence, some datasets encompass a great sampling coverage while others strive for a greater number of markers. In this study, we combined datasets to improve both sample sizes and number of markers to improve the power of analyses of tiger shark population structure. Comparing results from the datasets, they were not always congruent: the dataset with the greatest number of samples but the lowest number of molecular markers highlighted potential differentiation between the western Indian Ocean and the eastern Indian Ocean/Pacific, which was not identified with the dataset with larger number of markers but fewer samples. It demonstrates how difficult it may be to interpret results from multiple studies with different population samples and molecular screening effort.
This study investigated the population structure and demographics of the tiger shark by adding additional samples and new molecular markers to recent studies (Bernard et al., 2016;Holmes et al., 2017).
Compared to these previous studies, we carried out intensive sampling in the western Indian Ocean, a region in which only one location had previously been sampled (South Africa). We also used a higher number of microsatellite and mitochondrial loci in our analyses, thus expanding the picture of tiger shark population structure. The results confirmed here the genetic differentiation between the populations from the Indo-Pacific and the western Atlantic, with both microsatellite and mitochondrial markers, while an important genetic connectivity was detected within and between the Indian and the Pacific Oceans.
Furthermore, we investigated for the first time variations in effective population size at the scale of the Indo-Pacific and highlighted the probable occurrence of a bottleneck 2,000-3,000 years ago.

| Connectivity between the western Atlantic and the Indo-Pacific
Genetic differentiation between tiger sharks from the western Atlantic and the Indo-Pacific could only be investigated with the 8-msat and CR datasets as we did not obtain samples from the western Atlantic. By adding new sampling sites from the western Indian Ocean, the present study confirmed a genetic differentiation between the two regions. Both structure (using the LOCPRIOR model and at K = 3, i.e., not the first level of differentiation) and mitochondrial haplotype net- TA B L E 6 Characteristics of demographic parameter posterior distributions estimated using Galeocerdo cuvier individuals from Reunion Island with DiyaBc under Scenario 6 (population decrease at time t from N 1 to N b , followed by an expansion from N b to N 0 five generations ago) localities in the Indo-Pacific and the Atlantic compared to comparisons within the Indo-Pacific, both in the current study and in Bernard et al. (2016) and Holmes et al. (2017). Furthermore, the TCS haplotype networks constructed by Bernard et al. (2016) with the mitochondrial markers CR and COI clearly showed that no haplotypes were shared between the Indo-Pacific and the northwestern Atlantic, and only one was shared between the Indo-Pacific and Brazil, suggesting a real differentiation.

| Population connectivity within the Indo-Pacific
Patterns  Thus, altogether, weak genetic differentiation was highlighted between locations in the Indo-Pacific, both with microsatellite and mitochondrial markers. Tiger sharks are known to cross large oceanic expanses (Ferreira et al., 2015;Hammerschlag et al., 2012;Holmes et al., 2014;Werry et al., 2014), with records of transoceanic migrations in both the Indian  and the Atlantic Oceans (Afonso et al., 2017;Kohler et al., 1998;Kohler & Turner, 2001;Lea et al., 2015). These observations are thus in accordance with the weak genetic differentiation highlighted in the Indo-Pacific in this and previous studies.

| To be or not to be philopatric?
While the evidence remains scant, it is plausible that significant mitochondrial differentiation occurs. From their mitochondrial differentiation estimates, Bernard et al. (2016) inferred matrilineal population structure within the Indo-Pacific, which they linked to potential female site fidelity (i.e., philopatry) to reproductive areas, probably pupping sites. Female philopatry to nurseries has been hypothesized in many shark species, notably from discordances between mitochondrial and nuclear differentiation estimates (Karl, Castro, Lopez, Charvet, & Burgess, 2011;Pardini et al., 2001;Portnoy et al., 2014;Tillett, Meekan, Field, Thorburn, & Ovenden, 2012), but these discordances may also result from other processes (Prugnolle & de Meeus, 2002). Furthermore, the tiger shark seems to inhabit pelagic waters more frequently than other species such as the lemon shark Negaprion brevirostris, for which female philopatry to nurseries has been confirmed (Feldheim et al., 2013). The tiger shark may thus be less constrained in their movements and females may swim to various sites. It is noteworthy that no known tiger shark nurseries have been identified in the Indo-Pacific, leaving unanswered the hypothesis of female fidelity to coastal pupping areas. The tiger shark is also one of the few species for which multiple paternity (polyandry) has not been identified, although only eight litters in total have been investigated (Holmes et al., 2018;). Yet, this behavior is hypothesized to be linked to philopatry and more structured populations (Chapman, Prodohl, Gelsleichter, Manire, & Shivji, 2004), and the predominance of monoandry in the tiger shark may be another indication pointing to an absence of female philopatry to specific nurseries. All considered further studies, notably identification of nurseries and genotyping of females and juveniles for several years, as well as satellite tracking are needed to fully resolve the occurrence or absence of female tiger shark fidelity to nurseries.

| Low genetic diversity and bottleneck
The tiger shark displays moderate genetic diversity with a very low number of mitochondrial haplotypes and haplotype diversity for the sequences studied compared to other shark species, such as the bull shark Carcharhinus leucas , the great white shark Carcharodon carcharias (Pardini et al., 2001), the blue shark Prionace glauca (Veríssimo et al., 2017), the blacktip reef shark Carcharhinus melanopterus (Vignaud et al., 2014), or the tope shark Galeorhinus galeus (Chabot, 2015;Chabot & Allen, 2009). Furthermore, using the same protocol on the same date from extraction to marker testing for polymorphism, we characterized 20 microsatellite loci for the bull shark (Pirog et al., 2015), and only eight for the tiger shark (Pirog et al., 2016), which supports a lower genetic diversity in the latter species. Low genetic diversity has also been identified for another pelagic species, the basking shark Cetorhinus maximus (Hoelzel, Shivji, Magnussen, & Francis, 2006), for which only six haplotypes were identified using the CR marker, over its entire distribution range. This low genetic diversity was thought to be due to the occurrence of a bottleneck during the Holocene (Hoelzel et al., 2006).
Here, the Bayesian analysis, which combined both nuclear and mitochondrial information also provided evidence for the occurrence of a recent bottleneck experienced by the Indo-Pacific tiger shark population, 2,000-3,000 years ago (during the Holocene).
Approximate Bayesian computation has been proven to be one of the most accurate methods to detect recent demographic changes, as evidenced by the study of a large marine mammal for which historic fisheries data were available, the Antarctic fur seal Arctocephalus gazella (Hoffman, Grant, Forcada, & Phillips, 2011).
Estimates of effective population sizes for the tiger shark remained difficult to infer, but the ancestral size of this population was around 5,000 individuals (95% CI = [1,120,9,710]) and the bottleneck resulted in an effective population size as low as 111 individuals (95% CI = [43,369]). This bottleneck may well be responsible for the low genetic diversity presented by the species. Decreases in effective population size within the Holocene have also been identified for several large marine species, notably sea turtles (Molfetti et al., 2013), whales (Baker & Clapham, 2004), and elephant seals (de Bruyn et al., 2009), and also for terrestrial megafauna (Brook & Bowman, 2002). Notably, a similar signal to the one we identified here for the tiger shark was also found in the scalloped hammerhead shark Sphyrna lewini in its eastern Pacific range (Nance, Klimley, Galván-Magaña, Martínez-Ortíz, & Marko, 2011). Factors responsible for these population declines during the Holocene remain difficult to identify. This period is characterized by a general warming that likely induced population expansions of many marine species (Marko et al., 2010;Uthicke & Benzie, 2003), yet it is unclear how climate changes during this period may have led to population decreases of tropical and subtropical species. A widespread hypothesis for extinctions during the Pleistocene-Holocene is the emergence of diseases (Koutavas, Lynch-Stieglitz, Marchitto, & Sachs, 2002), notably for corals and sea urchins (Aronson & Precht, 2001;Lessios, 1988). Nevertheless, little is known about diseases and their impacts in sharks, despite an apparently robust immune system (Luer, Walsh, & Bodine, 2004;Walsh et al., 2006). One could also evoke the large volcanic eruptions that took place at this period. The ash layer with the increase of aerosols in the atmosphere absorbed solar radiation, leading to short-term cooling at regional to global scales, proportional to the magnitude of the eruptive episode (Robock, 2000;Sigl et al., 2015). The subsequent reduction in photosynthesis and the change in species distribution may have led to modifications of trophic networks, possibly impacting marine populations, notably apex predators. It is also possible that prehistoric fisheries were responsible for initiating this decrease, which was then intensified by modern fishing. Indeed, evidence of pelagic fishing was reported as early as 42,000 years ago in the western Pacific (O'Connor, Ono, & Clarkson, 2011). Other studies have provided evidence that prehistoric fishing may have impacted coastal ecosystems in parts of the world (Cooke, 1992;Erlandson & Rick, 2010). Furthermore, the presence of elasmobranchs (sharks, rays, and skates) in prehistoric fisheries is difficult to prove as their cartilaginous skeleton was rarely preserved in archeological sites (Rick, Erlandson, Glassow, & Moss, 2002). Thus, the impacts of prehistoric fishing on these populations may be underestimated. Regardless of the possible causes of tiger shark population decline in the Indo-Pacific, the extremely low effective population size estimated during the bottleneck (N b = 111) is lower than the minimal estimate of 500 thought to be required for a population to be able to adapt to environmental changes or pressures (Frankham, Briscoe, & Ballou, 2010). While this tiger shark population seems to have expanded since the bottleneck, we were not able to estimate the current effective population size, and this population may well be vulnerable in terms of genetic diversity. For now, the tiger shark is classified as Near Threatened in the IUCN Red List of threatened species (Simpfendorfer, 2009). Its conservation status might need to be revised according to the results highlighted in this and previous studies.

| CON CLUS ION
The study presented here confirms the genetic connectivity of tiger sharks within the western Indo-Pacific, using both microsatellites and mitochondrial markers. While mitochondrial differentiation estimates were slightly higher than those of microsatellite, further analyses are needed to confirm whether differentiation reflects female philopatry to specific nurseries. Individuals from the Indo-Pacific form a discrete population. Management and conservation programmes thus need to be designed at such scales in order to maximize potential efficacy. Meanwhile, localized intensive fishing may yet stand to impact the whole population.
Furthermore, the detection of low genetic diversity as well as a recent bottleneck (in the Holocene), during which the effective population size of tiger sharks in this region may have dropped as low as 111 individuals, points to a potentially resultant vulnerable population. Further assessments of the health status of this population as well as conservation plans are thus particularly needed to conserve the tiger shark within the western Indo-Pacific. Its conservation status might need to be revised toward a higher vulnerability level, as the ability of tiger sharks to withstand high levels of fishing pressure might be lower than previously thought.

ACK N OWLED G M ENTS
We acknowledge the Plateforme Gentyane of the Institut National

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

AUTH O R CO NTR I B UTI O N S
GC, EC, BJH, NEH, JEGN, AJT, PB, LV, and SJ provided samples. AP and HM did lab work. AP, HM, and VR analyzed the data. All authors contributed to the writing of the manuscript.HM and SJ designed research.