Introduction

The impacts of human activity on marine life are pervasive (O’Hara et al. 2021), and no ecosystem remains unaffected (Halpern et al. 2019). Human impacts extend far beyond coastal zones into the pelagic realm (Kroodsma et al. 2018) and consequently only 13.2% of the ocean can be considered refugia for wildlife (Jones et al. 2018). The cumulative effects of intensifying anthropogenic pressures such as climate change, pollution, habitat degradation and overexploitation threaten pelagic biodiversity, with overexploitation being the principal factor in the decline of oceanic sharks (White et al. 2019; Pacoureau et al. 2021). Oceanic shark abundance has depleted by 71% in the past half century, and they now represent one of the world’s most threatened groups of vertebrates (Dulvy et al. 2014; Pacoureau et al. 2021). The widespread and rapid decline of oceanic sharks is severely concerning due to their elevated risk of extinction (Field et al. 2009) and their critical functional roles in ecosystem connectivity, nutrient transport and fish community regulation (Ferretti et al. 2010; Heithaus et al. 2010; Pimiento et al. 2020).

The silky shark is a commercially important species of conservation concern (Cardeñosa et al. 2018; Hueter et al. 2018). Although one of the world’s most abundant sharks (Bonfil 2008), intense targeted and incidental fishing has substantially reduced their numbers by up to 54% over 45 years in some regions (Rigby et al. 2017). This reduction has been sufficient to warrant a raft of conservation measures, such as inclusion in Appendix 2 of the Convention on International Trade in Endangered Species of Wild Fauna and Flora (CITES) (2016) and uplisting to ‘vulnerable’ on the International Union for Conservation of Nature (IUCN) Red List of Threatened Species in 2017 (Rigby et al. 2017). Practical measures to reduce bycatch mortality, such as gear modifications and retention bans, have been adopted by multiple Regional Fisheries Management Organisations (RFMOs) (Serafy et al. 2012; Western and Central Pacific Fisheries Commission (WCPFC) 2013; Inter-American Tropical Tuna Commission (IATTC) 2019; Schaefer et al. 2019). Combined, however, these measures fall short to halt and reverse the ongoing global decline in silky sharks (Clarke et al. 2018). A growing focus for the conservation of mobile pelagic species such as silky sharks has been the designation of large marine protected areas (MPAs), in which extractive activities are prohibited or limited (Game et al. 2009; O'Leary et al. 2016; Dwyer et al. 2020). This growing need resonates with the global initiative to designate at least 30% of oceans as protected areas by 2030 (30 × 30) (Zhao et al. 2020). However, the applicability and effectiveness of static spatial protection to silky shark conservation remain unassessed and questionable due to their high mobility, vast home range and population structure (Le Quesne and Codling 2009; Mizrahi et al. 2018; Boerder et al. 2019; Kraft 2020). The effectiveness of MPAs for sharks is further limited by poor enforcement, inadequate governance and placement based on avoiding conflict with extractive or commercial uses rather than the actual ecological value of the site (Agardy et al. 2011; Edgar et al. 2014; Devillers et al. 2015; Jacoby et al. 2020).

The effectiveness of conservation measures for silky sharks is additionally constrained by limited knowledge of their distribution, abundance, and population structure, of which an improved understanding is urgently required as the foundation of successful temporal and spatial management strategies (Anderson and Jauharee 2009; Rigby et al. 2017; Hueter et al. 2018). The challenges of studying elusive species in the open ocean have led to a reliance on fishery-dependent data, which provide valuable information on geographic distribution but are limited by gear size selectivity, non-random distribution of fishing effort, and underreporting of shark landings (Kendall et al. 2009; Jacquet et al. 2010; Temple et al. 2019). The recent development of mid-water stereo-baited remote underwater video systems (BRUVS) has facilitated non-extractive sampling within previously inaccessible environments (Santana-Garcon et al. 2014; Letessier et al. 2019) and reduced our reliance on fishery-dependent data to enhance our understanding of pelagic wildlife (Letessier et al. 2017; Thompson et al. 2021).

We mined data from a mid-water BRUVS database curated by the Marine Futures Lab for 39 locations across three ocean basins to: (1) document the distribution, abundance and size structure of silky shark populations, (2) identify the most important environmental, geographic and anthropogenic predictors of their presence, abundance and body size, and (3) assess if their abundance is greater within designated MPAs than in locations not designated as MPAs. We focus our study on Exclusive Economic Zones (EEZs) to allow consideration of socioeconomic factors, which influence the distribution, abundance and diversity of marine wildlife (Brewer et al. 2012; Jaiteh et al. 2016) and the effectiveness of spatial protection (Mizrahi et al. 2018). Our results will allow insight as to whether MPAs have been designated in areas with remaining relatively high abundances of silky sharks, increase our ecological understanding of this vulnerable but poorly understood species and provide a baseline for future monitoring.

Methods

Survey methods

Our analysis is based on a dataset of video-based records of ocean wildlife collected using mid-water BRUVS. Each mid-water rig is composed of a lightweight 1.45-m long stainless-steel frame, a perpendicular 180-cm long bait arm, a 45-cm long perforated PVC canister containing 1 kg of pilchards (Sardinops sagax) and two high-definition small action cameras, in waterproof housings, that are mounted 80 cm apart and angled inwards by 8° on a horizontal base bar (Fig. 1). Stereo cameras enable identification, counts and length measurements of individuals present in the field of view within a range of ca. 10 m. Individual mid-water BRUVS are deployed in a longline configuration of ‘sets’ of three to five rigs, each separated by 200 m of horizontal line. Mid-water BRUVS are deployed, suspended at a depth of 10 m and then recovered following a 2-h soak time. All mid-water BRUVS were deployed during daylight hours to minimise the influence of crepuscular behaviour of fishes (Myers et al. 2016). Sampling design within locations was based on either stratified random sampling or generalised random tessellation stratification (Stevens and Olsen 2004) depending on survey purpose. All sampling followed standard protocols outlined in Bouchet et al. (2018) and was conducted under University of Western Australia ethics permit RA/3/100/1484 and appropriate Australian and international field permits.

Fig. 1
figure 1

Schematic of mid-water stereo-BRUVS in-water configuration. Adapted from Bouchet and Meeuwig (2015)

Video treatment and analysis

This analysis is based on video imagery from 1418 deployments of mid-water BRUVS during 66 surveys at 39 locations between April 2012 and March 2020 held in a dataset curated by the Marine Futures Lab at the University of Western Australia (Supplementary Information Appendix 1; Appendix 2). All locations were within the latitudinal range of silky sharks, between 42° N and 43° S (Ebert et al. 2013). Each longline of three to five mid-water BRUVS was considered and analysed as a single deployment given that the number of rigs per set varied and because the rigs on a given set are not statistically independent of each other (Hurlbert 1984). For each deployment, video analysts identified all observed individuals to the greatest taxonomic resolution possible, estimated the relative abundance of each taxon recorded as the maximum number of individuals in a given frame (MaxN), and measured individual fork lengths (FL). Each deployment was analysed for 120 min as a standardised sampling unit.

We first filtered the mid-water BRUVS dataset for all records of silky sharks. For consistency with the broader silky shark literature, in certain aspects of the analysis, we converted FL into total length (TL) with the equation: TL = −1.16 + 1.2 × FL (Branstetter 1987). We also allocated individuals into three life history classes: young-of-year (YOY), juveniles and adults. These categories were based on reports that silky sharks are born at sizes between 57 and 87 cm TL (Compagno 1984), and thus YOY were those individuals < 87 cm or less. To categorise adults, we assigned our sampling locations to each of the central Atlantic, central Indian, eastern Indian, eastern Pacific, western Pacific and southern Pacific Ocean regions, as sexual maturity is reached at different lengths depending on the varying life history characteristics of specific populations by region (Grant et al. 2018). Although there is a paucity of understanding surrounding the life history parameters of many pelagic shark, species and methods of estimation may influence observed differences in vital rates (Mukherji et al. 2021), and genetic analysis has revealed that silky sharks have multiple genetically distinct stocks within ocean basins (Kraft 2020). This genetic distinction likely contributes to intrinsic intraspecific variation in life history parameters, rather than these differences stemming from experimental design. Adults were classed as > 223 cm in the central Atlantic, > 180 cm in the eastern Pacific, > 194 cm in the western and southern Pacific and > 212 cm in the eastern and central Indian Ocean regions based on the review by Grant et al. (2018). Juveniles were classed as those individuals > 87 cm and less than the adult size of the region in which they were recorded.

Silky shark distribution and abundance

We evaluated the degree to which our mid-water BRUVS records reflect the predicted distribution of silky sharks from AquaMaps (www.aquamaps.org). AquaMaps is a tool for generating large-scale predicted range maps for marine species based on models which combine occurrence data from online databases with published literature (Kaschner et al. 2019). We first downloaded the probabilities of silky shark occurrence on a 0.5° × 0.5° grid based on depth (m), sea surface temperature (°C), salinity (psu), primary productivity (mgC·m−3·day−1) and distance from land (km) from AquaMaps. We then calculated the centroid of mid-water BRUVS sampling for each location within our dataset and associated mean relative abundance of silky sharks. Mean relative abundance was calculated as the average abundance of silky sharks per longline mid-water BRUVS deployment. The mean relative abundance for a location was then calculated as the average of the deployments in that location.

We visually compared the mean relative abundance of silky sharks with the predicted probability of occurrence from AquaMaps using QGIS 3.16 and generated a quadrant graph based on mean relative abundance and probability of occurrence. We also identified whether each location was within a designated MPA by referencing the MPAtlas (www.mpatlas.org). We used these outer boundaries of the MPAs rather than assigning locations to highly protected and partially protected zones within MPAs as, at the time of sampling, the majority of MPAs were too newly established to test zoning effectiveness. Moreover, as zoning can be strengthened or weakened, we felt that the outer MPA boundaries were the appropriate spatial unit for analysis.

To further understand our observed distribution and abundance of silky sharks, we tested for differences in abundance and total length among locations. Differences in abundance were tested using a Kruskal-Wallis test. Locations with sample sizes of fewer than five individuals were excluded to meet the requirements of the Kruskal-Wallis test (Riffenburgh 2006). Differences in total length were tested using a one-way ANOVA and Tukey’s post hoc test. Differences in abundance within MPAs compared to locations not designated as MPAs were tested using a non-parametric Mann-Whitney U test.

Drivers of silky shark presence, abundance and body size

We collated a set of 13 environmental, geographic and anthropogenic variables to be used as predictors of silky shark presence, abundance and body size (Supplementary Information Appendix 3). Variables were chosen based on known associations with marine wildlife and have been used previously in modelling predator distributions (Letessier et al. 2019). Environmental and geographic predictors were collated across a 0.01° × 0.01° grid. The environmental variables include mean sea surface temperature (SST) as a proxy for latitudinal patterns in species diversity universally observed across taxa (Tittensor et al. 2010) and its standard deviation (SSTSD) which is an indicator of frontal dynamics generating nutrient mixing and multilevel productivity (Edwards et al. 2013; Queiroz et al. 2016), and chlorophyll a (Chla) which indicates primary productivity and available trophic energy (Currie et al. 2004). These environmental variables were extracted from the NOAA Geophysical Fluid Dynamic Laboratory Earth System Model 2G (GFDL-ESM2G) and quantified as the annual mean for the sampling period. The geographic and geomorphological variables included depth as a dimension that fundamentally structures and constrains marine habitats vertically (Priede et al. 2013) and slope which provides an indication of habitat complexity (Friedman et al. 2012). We included three distance measures: (1) to coast as a measure of terrestrial energy availability (Gove et al. 2016) and a physical barrier restricting the horizontal extent of the marine habitat (Letessier et al. 2016); (2) to seamount, the presence of which is known to attract predators (Morato et al. 2010; Letessier et al. 2019); and (3) to reef as the distance to the centre of the Coral Triangle, the epicentre of fish diversity (Roberts et al. 2002).

Our anthropogenic variables included distance to human population, market and major port, as defined below. These distance metrics are proxies for cumulative effects of noise pollution (Nowacek et al. 2015), unreported fishing (Pauly and Zeller 2016), vessel strikes (Constantine et al. 2015), infrastructure development (Maire et al. 2016) and direct exploitation (Roff et al. 2018), all of which that can affect predators and have been previously employed as measures of human impact on sharks (D’agata et al. 2014; Mellin et al. 2016; Cinner et al. 2018; Juhel et al. 2019; Letessier et al. 2019). These three metrics also reflect some historical impacts that have occurred before the onset of modern record keeping (Pauly 1995; Pauly and Zeller 2016). Distance to the nearest human population was calculated using the LandScan 2016 database (Dobson et al. 2000). Distance to market was defined as the distance to the nearest human density centre using the World Cities spatial layer (ESRI) where density centres are defined as provincial capital cities, major population centres, landmark cities, national capitals and shipping ports. We included both distance metrics as small populations may impact predators locally (Bellwood et al. 2011) but impacts scale rapidly when industrialised markets are present (D’agata et al. 2014). The distance to the nearest medium or large port was calculated using the World Port Index database (ESRI) (Andrello et al. 2017).

We also included two socioeconomic indicators extracted at the level of sovereign state. The Human Development Index (HDI) captures elements of life expectancy, education and wealth. We used HDI 2015 values from the 2016 Human Development Report published by the United Nations of Development Program (UNDP). Fisheries dependency scores the importance of a coastal country’s fisheries in terms of their contributions to a national economy, employment and food security.

Prior to model fitting, we evaluated the multicollinearity between predictors using the Pearson coefficient of correlation and variance inflation factor (VIF). Except for one pair of predictors (log-distance to the coast and log-distance to the nearest population) and the standard deviation of sea surface temperature, all factors were weakly or moderately correlated (−0.55 > rs > 0.5) and had a variance inflation factor below five (Supplementary Information Appendix 4). We chose to keep all predictors except log-distance to the coast and standard deviation of sea surface temperature to result in 11 predictors out of 13.

We ran a first series of generalised linear models (GLMs) to predict silky shark occurrence, abundance and fork length using all predictors without interactions. To model abundance and fork length, only locations with shark occurrence were included. We employed presence-only models (Elith and Leathwick 2009; Gomes et al. 2018) to firstly manage the likelihood that absences may simply reflect the fact that sampling was conducted outside of the natural habitat range of silky sharks, and to secondly manage the zero-inflated nature of mid-water BRUVS occurrence records. We therefore investigated the drivers of silky shark presence, abundance and body size in areas in which they were observed to enhance our understanding of their current distribution. We used the binomial error distribution and a logit link function to predict the probability of silky shark presence, and we then used the Gaussian distribution to predict the log-abundance of sharks and their body length. We subsequently estimated the relative importance of each predictor in each model by summing Akaike weight values across all models that include the predictor. These summed Akaike weights (WAIC) range from 0 to 1, hence providing a means for ranking the predictors in terms of information content. We also estimated the significance of each predictor (F-value) in each model. We then ran similar GLM models but including all two-term interactions between pairs of predictors to assess whether predictions can be improved with higher model complexity. Models with and without interactions were compared using AIC values and Likelihood Ratio Tests (LRT). The performances of the binomial GLM models, predicting shark occurrence, were assessed using the area under the curve (AUC) statistic, for which values are considered random when they do not differ from 0.5, poor when they are in the range 0.5–0.7, and useful in the range 0.7–0.9, and excellent above 0.9. For Gaussian GLM models, we used the R2 values between observed and modelled values. We then tested for spatial autocorrelation in model residuals using Moran’s I test. Finally, we built partial regression plots showing the effect of each predictor while controlling for all the others on the probability of detecting silky shark, their mean log-abundance and their mean fork length. All analyses were performed in R version 3.5.1 using the ‘glm’, ‘dredge’ and ‘visreg’ functions.

Results

Silky shark distribution

The dataset contained 945 silky sharks on 413 deployments from 18 locations (Table 1). An examination of AquaMaps with our samples superimposed shows that no silky sharks were observed where they were not expected to be (Fig. 2). However, abundance was not correlated with probability of occurrence (Fig. 2). High abundances of silky sharks only occurred at three locations (8%) with a high probability of occurrence (probability > 0.5). Conversely, low or zero abundances were observed at 10 locations (26%) with a low probability of occurrence (probability < 0.5). Interestingly, 26 locations (66%) were categorised as having a high probability of occurrence yet were characterised by low or zero abundance of silky sharks. Of the 26 high probability-low abundance locations, 17 were within MPAs.

Table 1 Location, including ocean region, of all silky shark (Carcharhinus falciformis) records, whether the location is within a marine protected area (MPA), the year of effective MPA enforcement extracted from the MPAtlas, the number of BRUVS deployments, total number of individuals observed and measured, the mean total length (TL) and associated standard error (SE), and minimum and maximum TL. The * indicates that sampling occurred prior to implementation of spatial protection
Fig. 2
figure 2

Observed distribution of silky shark (Carcharhinus falciformis) relative to predicted abundance from AquaMaps represented (a) spatially and (b) as a quadrant plot function of probability relative to observed relative abundance by location. Locations were plotted against four categories: (1) low probability of occurrence and high abundance, where low probability of occurrence is < 0.5; (2) low probability of occurrence and low abundance; (3) high probability of occurrence and low abundance, where high probability of occurrence is > 0.5; (4) high probability of occurrence and high abundance. Sampling effort is indicated by diameter of location marker and represents the number of mid-water stereo-BRUVS deployments

Silky shark abundance

Silky shark abundance varied significantly across locations (H (11) = 131.08, p < 0.001, n = 934) with a mean relative abundance of 0.26 ± 0.03 SE MaxN per mid-water BRUVS deployment (0.01–1.90 MaxN per mid-water BRUVS deployment, n = 945). We observed the greatest mean relative abundances at Clipperton Island (1.91 ± 0.23 SE MaxN per mid-water BRUVS deployment), Malpelo Island (1.59 ± 0.59 SE MaxN per mid-water BRUVS deployment) and Revillagigedo Islands (1.41 ± 0.17 SE MaxN per mid-water BRUVS deployment), which are all located in the eastern Pacific Ocean. In contrast, we observed the lowest mean relative abundances at Perth Canyon (0.003 ± 0.002 MaxN per mid-water BRUVS deployment), Shark Bay (0.003 ± 0.002 MaxN per mid-water BRUVS deployment) and New Caledonia (0.007 ± 0.006 MaxN per mid-water BRUVS deployment), located in the eastern Indian and western Pacific Oceans, respectively.

At the time of sampling, 20 of the 39 locations were designated as MPAs (Supplementary Information Appendix 1). MPAs had a median age of 6.5 years. The sites at which sampling took place around Ascension Island, the Azores, the Cocos (Keeling) Islands, Niue, Palau, the Selvagens Islands and Tristan da Cunha gained spatial protection after the sampling period (Supplementary Information Appendix 1).There was no significant difference in the abundance of silky sharks between MPAs (n = 20) and locations not designated as MPAs (n = 19) (W = 226, p = 0.28).

Silky shark body size

We estimated FL for 511 (44%) individuals with the remainder not measured due to distance from or poor orientation of the individual shark to the BRUVS. The mean FL was 141 cm ± 0.07 SE (range of 57–275 cm). The mean converted TL was 168 cm ± 0.09 SE (range of 68–329 cm). Most individuals ranged in size between 140 and 150 cm TL (Fig. 3). Mean TL of silky sharks varied significantly among locations (ANOVA, F (11, 490) = 17.42, p < 0.001, n = 502) and among regions (ANOVA, F (4, 497) = 26.17, p < 0.001, n = 509). Silky sharks in the eastern Indian Ocean region had a significantly larger mean TL than individuals in the central Indian Ocean region (p < 0.01), and those in the eastern Pacific Ocean region had a significantly larger mean TL than individuals in the western Pacific Ocean Region (p < 0.001) (Fig. 3).

Fig. 3
figure 3

a Length frequency distribution of silky sharks (Carcharhinus falciformis) (n = 513) and b variation in mean total length (TL) across ocean regions where they were recorded, with error bars displaying standard error (n = 512). The southern Pacific region was omitted from the graph due to only being able to measure one individual

Juveniles were the most frequently observed life stage with 377 records (74%), while 125 individuals were classified as adults (24%), and 9 individuals were classified as YOY (2%). The YOY were recorded at four locations; juveniles were recorded at 14 locations; and adults were recorded at 17 locations (Fig. 4). There was no strict spatial segregation by life history stage since YOY co-occurred with juveniles and adults at 10% of locations. Juveniles were exclusively present at one location (Montebello Islands in the eastern Indian Ocean) and adults were exclusively present at four locations (New Caledonia in the western Pacific Ocean, Niue in the southern Pacific Ocean, and Perth Canyon and Shark Bay in the eastern Indian Ocean). There were no locations where YOY were exclusively present.

Fig. 4
figure 4

Locations where young-of-year (YOY), juvenile and adult silky sharks (Carcharhinus falciformis) were recorded by stereo-BRUVS. Sampling effort is indicated by diameter of location marker and represents the number of mid-water stereo-BRUVS deployments

Drivers of silky shark presence, abundance and body size

Our GLM models accurately predicted the presence, abundance and body size of silky sharks (Table 2), and in the cases of abundance and body size, included at least one measure of human pressure among the most important factors. Models with interactions performed significantly better than no-interaction models for shark abundance and body size but it was the opposite for shark presence (Table 3). A full report of the parameters and their interactions is available in Supplementary Information Appendix 5. For all models, the Moran’s I statistic was not significant (p > 0.05) showing that there is no spatial autocorrelation in the residuals. To highlight the main effects of predictors, we only present no-interaction model outputs.

Table 2 Results from generalised linear models (GLMs) predicting the presence, the mean log-abundance and mean fork length of silky sharks (Carcharhinus falciformis) as a function of two environmental, four geographic and five human factors
Table 3 Comparison of generalised linear models (GLMs) predicting the presence, the mean log-abundance and mean fork length of silky sharks (Carcharhinus falciformis) without any interaction between the 12 factors or with all second order interactions

The presence of silky sharks was best predicted by the model with an AUC value of 0.82 (Table 3). Variation in presence was primarily explained by SST, depth, distance to reef, distance to seamount and chlorophyll a (Fig. 5). The likelihood of presence increased with warmer SST, decreasing depth and greater distance to reef, and decreased with greater distance to seamount and increasing chlorophyll a.

Fig. 5
figure 5

Silky shark (Carcharhinus falciformis) presence as a function of environmental, geographic and human factors predicted using generalised linear models

The model for silky shark abundance explained 53% of variation (Table 3). Variation in abundance was primarily explained by distance to port, depth, fisheries dependency, HDI, and distance to seamount (Fig. 6). Abundance increased with greater distance to port and decreasing depth, and decreased with higher fisheries dependency, higher HDI and greater distance to seamount.

Fig. 6
figure 6

Silky shark (Carcharhinus falciformis) log mean abundance as a function of environmental, geographic and human factors predicted using generalised linear models

The model for silky shark fork length explained only 29% of variation with no interactions but 66% with interactions (Table 3). Variation in body size was primarily explained by distance to port, HDI, distance to population, distance to seamount, SST and distance to reef (Fig. 7). Body size increased with greater distance to port, higher HDI, greater distance to population and greater distance to seamount, and decreased with warmer SST and greater distance to reef.

Fig. 7
figure 7

Silky shark (Carcharhinus falciformis) mean fork length as a function of environmental, geographic and human factors predicted using generalised linear models

Discussion

This rare opportunity to study a cosmopolitan mobile predator on a large spatial scale yields new insights into their distribution, abundance and demographics. BRUVS have facilitated marine management through enabling non-extractive sampling within previously difficult-to-access environments (Letessier et al. 2014; Santana-Garcon et al. 2014; Letessier et al. 2019). The versatility of BRUVS in their application to both seabed and mid-water environments is demonstrated at the species and ecosystem-level, with insights gained into the impact of anthropogenic pressures on shark abundance (Goetze et al. 2018), size (Letessier et al. 2019) and behaviour (Juhel et al. 2019). Our study represents the first application of mid-water BRUVS to the silky shark across such a broad spatial scale and the creation of a baseline mid-water BRUVS dataset provides a valuable foundation for ongoing monitoring and management.

The presence, abundance and body size of silky sharks can be predicted by a parsimonious set of environmental, geographic and anthropogenic factors. The outcomes of our models were mainly expected and consistent with previous findings (Lezama-Ochoa et al. 2016; Letessier et al. 2019; Lopez et al. 2020; Thompson et al. 2021). Temperature and depth were the primary environmental and geographic factors driving presence and abundance, with silky sharks being most abundant in the relatively warm and shallow waters of the epipelagic zone. Indeed, previous studies have revealed that silky sharks spend nearly all of their time in the top 200 m of the water column in temperatures of 24–30 °C (Hueter et al. 2018; Hutchinson et al. 2019; Curnick et al. 2020). Temperature was also an important factor in predicting body size, with smaller individuals associated with warmer temperatures. The effect of temperature on the size of marine ectotherms is well-documented (Bergmann 1847; Audzijonyte et al. 2019), with the Gill Oxygen Limitation Theory arguing that smaller individuals with lower oxygen demands are able to survive in warmer, less oxygenated waters (Pauly 2019). Our study consequently demonstrates the applicability of environmental preference modelling, based on physiology, in identifying areas of high silky shark probability at different life stages.

Our models also demonstrate the impacts of humans on silky shark abundance and body size. The contribution of human footprint in our models is unsurprising given the intensified fishing pressure and habitat degradation associated with human settlements (Ward-Paige et al. 2010; Cinner et al. 2018). Proximity to port was an important factor in determining abundance, mirroring the findings of previous studies that document lower abundances of sharks in less remote areas (Goetze et al. 2018; Juhel et al. 2019; Letessier et al. 2019; MacNeil et al. 2020). Our documented size structure of observed silky sharks can additionally provide insight into the influence of human pressures on this trait. Variation in mean length across ocean regions may reflect differences in the life history parameters of distinct populations resulting from natural variation (Grant et al. 2020). Alternatively, variation in body size may reflect heterogeneity in fishing pressure across locations as overfishing typically leads to populations with reduced size (Tu et al. 2018). As our models identified proximity to population and port as primary drivers of silky shark body size independent of ocean basin, fishing pressure is likely to affect silky shark size and consequently population life stage structure. This is consistent with data for other oceanic shark species with documented declines in the average length of great hammerheads (Sphyrna mokarran), scalloped hammerheads (Sphyrna lewinii) and tiger sharks (Galeocerda cuvier) by 16–22% from 1997 to 2017 (Roff et al. 2018) and of whale sharks (Rhincodon typus) by nearly 2 m from 1995 to 2004 (Bradshaw et al. 2008). Given the importance of age-structured demographic models in the management of silky sharks (Grant et al. 2020), further investigation is required to establish whether our results reflect variation in life history characteristics, seasonality, fishing effort or study design.

We did not observe the expected correlation between probability of occurrence and abundance. Reassuringly, we did not observe silky sharks where they had a low probability of occurrence which validates the underlying AquaMap. Our observation of no or low silky abundances in high probability locations may reflect a number of alternatives. The large cluster of high probability-low abundance locations may simply indicate regional heterogeneity of silky shark density, with individuals aggregating outside of our sampling sites. Our findings may also reflect a limitation of our sampling as mid-water BRUVS are suspended at 10 m, which may reduce their likelihood of detection on BRUVS. The common depth of silky sharks is reported between 0 and 500 m (Last and Stevens 1994) and while they are known to frequent <50 m (Hueter et al. 2018; Hutchinson et al. 2019), they have been recorded at a maximum depth of 1112 m (Curnick et al. 2020). Finally, our observed abundances could reflect real declines in silky shark numbers in these locations. Without baseline BRUVS data we are unable to conclude that abundance has decreased across our sample locations. However, our findings that silky sharks appear absent from where they should be, have low abundance even in areas of highly suitable habitat, and have juvenile-dominated size structures, all suggest that silky shark populations are heavily depleted. We also demonstrate that presence, abundance and body size respond strongly to human impacts. Our findings corroborate the wider literature which reports ongoing declines in abundances of silky sharks and other pelagic sharks (Pacoureau et al. 2021). Within the past century there have been estimated regional declines between 50 and 79% in shortfin mako (Isurus oxyrinchus) abundance (Rigby et al. 2019a) and over 80% in great hammerhead and smooth hammerhead (Sphyrna zygaena) shark abundances (Rigby et al. 2019b; Rigby et al. 2019c), primarily driven by fishing.

Our findings reinstate the importance of identifying and protecting critical habitats for silky shark conservation. Proximity to seamounts was a fundamental driver of presence and abundance, mirroring previous documentation of the strong association between silky sharks and seamounts (Melo et al. 2003; Morato et al. 2010; Letessier et al. 2014; Thompson et al. 2021). Interestingly, we found that proximity to seamounts was associated with smaller body size, perhaps suggesting that seamounts provide refuges for juvenile silky sharks (Orr 2019) and further demonstrating the importance of these topographical features as key habitats for silky sharks (Letessier et al. 2019). Seamounts are thought to aggregate pelagic species through a combination of their increased productivity and magnetic signatures which provide navigational aids for migratory species (Holland and Grubbs 2007; Litvinov 2007). Tagged silky sharks associated with a seamount complex in the British Indian Ocean Territory, one of the sampling locations within this study, exhibited diurnal movement patterns speculated to be related to feeding (Curnick et al. 2020). Spatial protection of seamounts may represent an important avenue for the conservation of silky sharks and other migratory species as they serve as a landmark in a homogeneous pelagic environment and their clear boundaries enable more practical surveying and enforcement (Morato et al. 2010; Jacoby et al. 2020).

Nursery sites also represent important habitats for sharks due to the protection from predation that they offer to pups (Vandeperre et al. 2014). Protection of nurseries is consequently an important component of silky shark conservation (Kinney and Simpfendorfer 2009) as juvenile survivorship underpins the ability of overexploited populations to recover (Heithaus et al. 2008). Our observed distribution of YOY provides insight into the nature of silky shark nurseries. Nurseries are characterised by spatial segregation from adults (Heupel et al. 2019) which we did not observe. This suggests that either we did not observe any ‘true’ nurseries or that this generality does not apply to silky sharks. It is also possible that while our overall database is large, the number of YOY is too few to make inferences about nurseries. Despite the lack of life history segregation, the locations where we observed YOY did share some commonalities that enhance pup survival. Shark nurseries are typically characterised by warm waters, heterogenous bathymetry, high productivity and relatively high abundance of potential prey species (Heupel et al. 2019) as is the case at our Malpelo Island, Palau, Wandoo and British Indian Ocean Territory locations. Further research into the fine-scale environmental preferences of silky shark pups would facilitate the development of models specifically identifying potential nursery sites similar to that conducted by Oh et al. (2017).

Finally, we recorded markedly low abundances of silky sharks within established MPAs, and there was no observed difference in abundance between locations within MPAs and locations where no MPA had been established. This suggests two alternatives. First, given the median ‘life’ of these MPAs of 6.5 years prior to sampling, the effectiveness of MPAs in conferring protection to silky sharks may be masked by this species’ slow maturation and low fecundity, resulting in a time lag between MPA enforcement and increased abundance. The age of MPAs (>10 years) has been identified as one of the five key factors in determining their success (Edgar et al. 2014). Thus, the benefits of the small increase in global MPA coverage from 1.9 to 8.2% since 2010 (O’Leary et al. 2018; Marine Conservation Institute 2023) in response to international spatial protection targets is unlikely to be effective at short term. Second, it may be that MPAs have simply not been placed in areas where silky sharks are abundant. MPAs are often planned with limited knowledge of the spatial ecology of highly mobile species and are consequently unlikely to encompass critical habitats (Daly et al. 2018). It is reassuring that Ascension Island’s entire EEZ was designated as an MPA in 2019 following our sampling, given we observed over one-third (38%) of all individuals there. Our study identifies the Eastern Tropical Pacific Ocean as a hotspot for silky sharks, with the greatest relative abundances of silky sharks observed at Clipperton Island, Malpelo Island and the Revillagigedo Islands. Again, it is reassuring that Clipperton Island gained spatial protection in the year of sampling and the marine reserve surrounding the Revillagigedo Islands was expanded the year following sampling, and we recommend ongoing monitoring to evaluate the effect of MPA establishment and expansion on silky shark abundance. Our findings additionally support the continued management and strengthening of the Eastern Tropical Pacific Marine Corridor (CMAR), which encompasses Malpelo Island and the Galápagos Islands (Enright et al. 2021), to ensure adequate protection of silky sharks transiting between these locations.

The first caveat of our study relates to the detection of silky sharks. We may have overestimated the abundance of silky sharks if individuals are recorded on multiple mid-water BRUVS within a single longline deployment. This problem can be addressed by the use of MaxInd, a metric of the maximum number of animals where individuals can be reliably identified by markings (Sherman et al. 2018). However, Thompson (2021) applied MaxInd to blue sharks (Prionace glauca) and found that the results of MaxInd and MaxN were indistinguishable in terms of overall patterns. We thus suggest that overall, the use of MaxN, averaged across the longline deployment, is informative in terms of silky shark distributions. Conversely, we may have underestimated the abundance of silky sharks owing to the aforementioned inherent limitation of deploying mid-water BRUVS at 10 m. A potential solution to determine whether observed absences were true absences or simply an inability of mid-water BRUVS to detect silky sharks would be to cross-reference mid-water BRUVS data with fisheries-dependent data to confirm if silky sharks are present beyond the 10 m depth limit. However, employing fisheries-dependent data is challenging due to their poor taxonomic resolution, particularly within the Carcharhinus genus (Cashion et al. 2019; Henderson 2020). A natural progression of this work is the validation of mid-water BRUVS data with fisheries-dependent data, which is especially pertinent given the expansion of monitoring pelagic wildlife with mid-water BRUVS as part of the UK government’s Global Ocean Wildlife Analysis Network.

The second caveat relates to the scope of our sampling: we recognise that no mid-water BRUVS were deployed in developing countries which are within the range of silky sharks and we are missing key locations such as the African continent. Our sampling was also limited to EEZs and does not extend to the high seas. Silky sharks are circumglobally distributed in pelagic waters and thus knowledge of their distribution and abundance in the high seas would be informative. We also recognise that silky sharks are highly mobile and therefore transit between EEZs and across the high seas, and between MPAs and unprotected areas. This has two implications; firstly, it remains unclear whether even large MPAs within EEZs can recover silky shark populations (Game et al. 2009; Le Quesne and Codling 2009; Moffitt et al. 2009). Secondly, it is necessary to consider what management measures apply to silky sharks outside of protected areas. Silky sharks are classed as no-retention species by some of the Regional Fishery Management Organizations (RFMOs) operating across several of our sampling locations (see Cardeñosa et al. 2020 for a review). However, the initial trauma of capture can still cause mortality (Tolotti et al. 2015), with post-release mortality of silky sharks incidentally captured in purse seine fisheries reported as up to 84% (Hutchinson et al. 2015) and post-release mortality after longline capture reported as 20–38% (Musyl and Gilman 2018; Musyl and Gilman 2019). Moreover, detailed statistics on post-release mortality and vessel compliance are poor (Clarke et al. 2013; Musyl and Gilman 2019) and thus incorporating the influence of retention bans is problematic. A comprehensive review of the effectiveness of all conservation measures which apply to silky sharks would be a valuable avenue of future research.

Our novel application of mid-water BRUVS to assess silky shark presence and abundance exemplifies the value of this tool in increasing our understanding of pelagic species inhabiting previously difficult-to-access environments. Whilst environmental and geographic characteristics were identified as important determinants of silky shark presence, abundance and body size, the strong response of these attributes to human pressures exemplifies the pervasive impact of anthropogenic activities on large pelagic predators and supports an expansion of the current MPA network to encompass critical habitats such as seamounts. We provide a platform for ongoing monitoring of silky shark populations to investigate the applicability of static spatial protection to silky sharks and other highly mobile species.