Evaluating the use of marine protected areas by endangered species: A habitat selection approach

1. Optimizing the design of marine protected area (MPA) networks for the conservation of migratory marine species and their habitats involves a suite of important considerations, such as appropriate scale requirements and the distribution of anthropogenic impacts. Often, a fundamental component of the conservation planning process is delineating areas of high use or high biodiversity within a region of interest. 2. However, basing conservation strategies off merely the number of individuals in an ecosystem is outdated and potentially subject to arbitrary thresholds. To be effective at protecting marine megafauna, MPAs would ideally encompass habitats used by focal species. Through satellite-tracking studies, evidence of whether species actually use protected areas is emerging. 3. Here, we present a multispecies perspective on habitat selection within existing MPAs throughout the Floridian ecoregion, which encompasses coastal Florida and the Gulf of Mexico. Using an 11-year satellite-tracking dataset on 235 marine turtles, we used integrated step selection analysis to quantify the effects of sea turtle behavioural state (identified by a switching state-space model), protected area status, chlorophyll and bathymetry on habitat selection. 4. Our results show that sea turtles do select for existing protected areas, specifically multi-use zones, while controlling for the effects of depth and primary productivity. However, our analysis revealed that turtles showed no selection for the no-take zones within MPAs, during either transiting or foraging. 5.ThesefindingscontributetotheexistingliteraturebaseofMPAuseforhighlymobile, imperilled species and could inform management of existing MPAs or changes to zoning, specifically multi-use to no-take. Our use of a robust spatial modelling framework to evaluate habitat selection relative to MPAs could be incorporated into conservation planning to build MPA networks designed to accommodate migratory species.


INTRODUCTION
Marine environments have been subjected to an increasingly diverse set of human impacts (i.e. commercial and recreational fishing, infrastructure development, pollution, shipping and mining explorations) that have contributed to a steady decline in ocean ecosystem health (Halpern, Selkoe, Micheli, & Kappel, 2007;Klein et al., 2015). In response to these impacts, marine protected areas (MPAs) have become a widespread conservation strategy employed to protect and manage marine biodiversity as well as support sustainable use of ocean resources (Roberts, Valkan, & Cook, 2018). Internationally recognized conservation targets, such as the Convention on Biological Diversity (CBD) Aichi Target 11 and the Sustainable Development Goal 14.5, illustrate the importance of MPAs on a global scale (Campbell & Gray, 2019). Furthermore, these targets emphasize a focus beyond simple area-based metrics to ensure MPAs are designed to be ecologically connected, representative of vulnerable habitats and ecosystems, and effectively managed (Roberts, Duffy, & Cook, 2019; Secretariat of the CBD, 2011).
MPAs can be evaluated based on numerous indicators of performance relative to either management, economic or biological objectives. Methods of assessing MPA effectiveness for biodiversity conservation include evaluating coverage of representative habitats (Roberts et al., 2019), consistent monitoring of species richness or abundance (Strain et al., 2019), quantifying the benefit to fisheries due to a spillover effect (Di Lorenzo, Claudet, & Guidetti, 2016) or defining hotspots within or surrounding current boundaries to focus future conservation efforts (Worm, Lotze, & Myers, 2003). As MPAs are established under a broad range of contexts, results and recommendations from MPA evaluation studies are highly variable. To ensure effective design and placement of future MPAs as well as favourable outcomes for vulnerable biodiversity, it is important to understand the latest advances in impact evaluation methodologies (Ahmadia et al., 2015).
The successful management and conservation of endangered species relies on an accurate understanding of their distribution, movement patterns and interactions with their environment (Jeffers & Godley, 2016). Advances in satellite tracking of individuals in both marine and terrestrial environments have paved the way for improved ecosystem-based management approaches, capable of integrating food web interactions, life history processes and environmental parameters into conservation planning (Hays et al., 2019;Trathan et al., 2018). Telemetry data can reveal important habitats for key species across wide regions, as many marine species undertake long distance migrations throughout their life spans .
In order to effectively study this movement behaviour, models have been designed to account for multiple levels of biological and statistical complexities, such as irregular time intervals and location error, often associated with satellite telemetry data (Auger-Methe et al., 2016;Hoover et al., 2019;Jonsen, Flemming, & Myers, 2005).
Often, satellite tracking efforts can reveal high-use areas for foraging or breeding which can be summarized for changes in conservation policy (Hays et al., 2019;Lea, Humphries, Brandis, Clarke, & Sims, 2016). Delineating areas of high use or high biodiversity, commonly referred to as hotspot analysis, has become a fundamental component of conservation planning due to its feasibility and cost-effectiveness, as protecting the full range of biodiversity is never a realistic target (Marchese, 2015;Sussman et al., 2019). While typically designated on a case-by-case basis, hotspots are broadly defined as geographic areas with persistent high levels of either species abundance, richness or endemism (Possingham & Wilson, 2005). The spatial scale of hotspots, both in terms of geographic extent and resolution, can vary substantially depending on the sampling methodology and the technology being utilized (Hazen et al., 2013;Possingham & Wilson, 2005).
Hotspots in marine systems can be defined using qualitative or quantitative criteria. Qualitative approaches, such as mapping abundance and comparing densities visually, are easy to adapt regardless of species or ecosystem but are limited because often they do not reflect temporal change (Marchese, 2015). Additionally, qualitatively defined hotspots can be subjective and based on arbitrary thresholds that do not accurately reflect long-term data or actual use of the specified area by target species (Piacenza et al., 2015). Quantitative approaches to detect hotspots, such as kernel density estimation and the Getis-Ord G i * statistic, are spatially explicit and generally provide more consistency to inform long-term conservation agendas (Harvey, Nelson, Fox, & Paquet, 2017). However, significant variation can exist among statistical methods which can lead to misidentifying some areas or vastly different results across models, depending on the input parameters.
The hotspot approach is frequently implemented globally in management or conservation regimes as a result of telemetry studies (Ceballos & Ehrlich, 2006). The dynamic nature of the marine environment, due to complex physical processes, means that boundaries and features are constantly shifting. As a result, the effectiveness of conservation strategies reinforced by telemetry-based marine hotspots is often the subject of debate (Marchese, 2015). Here, we propose an alternative, spatially explicit method to the hotspot approach, designed to evaluate the use of MPAs by three species of endangered sea turtle: loggerhead (Caretta caretta), green (Chelonia mydas), and hawksbill (Eretmochelys imbricata). We used integrated step selection analysis (iSSA) to measure habitat selection, with protected status as one of the habitat covariates. Integrated SSA jointly estimates the relative probability of selection of different habitat covariates as well as the selection-free movement properties of the animal's trajectory (Avgar, Potts, Lewis, & Boyce, 2016). SSA is a special case of resource selection analysis (RSA), which use a use versus availability design to estimate the strength of habitat selection. While the step selection framework has previously been utilized in terrestrial studies (Abrahms et al., 2016;Signer, Fieberg, & Avgar, 2019;Thurfjell, Ciuti, & Boyce, 2014), it is rare that this tool is applied in marine research.

Study area
We constrained our analysis to the Floridian ecoregion, defined by which are subject to change depending on the jurisdiction. All MPAs included in this analysis were divided into two categories, no-take and multi-use, based on the specifications outlined in the zoning scheme or management plans. MPAs not defined by any zoning scheme were classified as multi-use in our analysis, unless specifically designated entirely as no-take.

Data collection
Throughout the 11-year tracking effort, sea turtles were tagged after capture either in-water or after a nesting event in the Dry Tortugas National Park (Florida), Everglades National Park (Florida), Broward County (Florida), or Gulf Shores (Alabama) (Supporting Information, Figure SA1). Individuals were captured and tagged following identical methods to those outlined in previous research (Hart, Lamont, Fujisaki, Tucker, & Carthy, 2012). Platform transmitter terminals were adhered to the turtle carapace using slow-curing epoxy (two-part Superbond epoxy). These tags were programmed to collect location data 24 hr per day and transmit once a day through the Argos satellite system when the turtle breached the surface of the water. Beginning in 2011, tags on nesting loggerheads were programmed to transmit every third day from November 1 through April 1 in order to preserve battery life. Raw location data were downloaded from the Wildlife Computers portal.

Data analysis
Wildlife movement data collected via marine satellite tags typically have two properties that make them difficult to analyse: (1)  The SSM fills gaps in the existing data by estimating a correlated random walk based on observed locations. When gaps of many days occur in the raw data, the resulting track becomes less informed by data and thus less reliable the longer the gap is. Thus, we split individual tracks with gaps longer than 25 days, and we passed these tracks to the SSM as if they were separate individuals. After fitting the SSM, we recombined the modelled daily locations for each turtle. The SSM categorized each step into one of two behavioural categories: one which was characterized by relatively shorter step lengths and sharper turning angles, and another which was characterized by relatively longer step lengths and straighter turning angles. We interpreted the behaviour corresponding to the former as 'area-restricted search' (ARS) and the latter as 'transiting' . The SSM output was then used as the input for the iSSA.
For the iSSA, observed steps are compared to available steps which are randomly sampled from proposed theoretical distributions of step lengths and turning angles (Avgar et al., 2016). We chose to sample step lengths from a gamma distribution and the turning angles from a von Mises distribution. Because the available locations are considered at the step scale, the iSSA framework is ideal for estimating selection of MPAs in the face of tagging bias (i.e. when animals are tagged within or near the MPA that we wish to evaluate). Specifically, animals tagged within the MPA will likely have a high proportion of their used locations within that MPA, simply due to autocorrelation. Any analysis that compares use to availability (e.g. traditional habitat selection analysis) without controlling for this dependency would conclude the animal selects for the MPA. But in the iSSA framework, availability is defined at the step level (rather than the study area or home range). Those animals tagged within the MPA will also have a high proportion of available locations within that MPA, thus eliminating any bias in the estimation of habitat selection. Because iSSA can take the movement process into account, it also has the capacity to account for autocorrelation in the data and results in habitat selection inference unbiased by the movement process (Avgar et al., 2016).
We generated 15 random steps for each observed step and analysed the resulting dataset of used and available steps using case-control logistic regression (also called conditional logistic regression), where each observed step is matched with its corresponding random steps (Duchesne, Fortin, & Courbin, 2010;Fortin et al., 2005). We fit one model for all the ARS steps and one model for all the transiting steps, thus parameterizing two integrated step selection functions (iSSFs).
We ran the iSSA using the function 'fit_issf()' from the R package 'amt' (Signer et al., 2019). Each iSSF included the same six covariates: (1)  We assigned the covariates to each step based on their value at the endpoint of the step. We assigned protected status to each step using a shapefile from the National Marine Protected Areas Center (https://marineprotectedareas.noaa.gov/). We assigned depth by using the R package 'marmap' to query the data from NOAA and attach it based on location (Pante & Simon-Bouhet, 2013 Figure 3). A positive log-RSS value means the animal selects for x1 relative to x2, a negative log-RSS value means the animal selects x2 relative to x1 and a log-RSS of 0 means it does not select either over the other.
Finally, we generated predicted relative selection maps in ArcMap to visualize the utilization distribution where selection is conditional relative to fixed covariates (Supporting Information, Figure SA2). We designated selection to be relative to unprotected waters in average (18 m) water depth. Additionally, as it was difficult to account for temporal variability at this scale, chlorophyll was held constant at its mean.

State space modelling
After filtering extreme outlier locations, we were left with 303,200 Argos locations. We split 36 individual tracks at gaps of greater than 25 days, and some tracks were split more than once due to multiple large gaps. One individual (tag ID 175681), tagged in the Dry Tortugas, was hindering the fit of the SSM and was removed from further analysis. The final number of tracks used as separate IDs in the SSM was 275.
After fitting the SSM, we were left with a total of 44,329 daily locations. Individual turtles had a mean of 161 daily locations (SD = 135, range = 11-1170; Supporting Information, Figure SA3). The majority of locations (95.6%) were assigned to the behavioural state 'ARS' .

Step selection analysis
We converted the locations from the SSM to steps, generated the ran-  (Figure 3a). They avoided areas of the highest chlorophyll concentrations while in ARS but showed no selection for chlorophyll while transiting (Figure 3b). A detailed summary of step selection coefficients is provided in Supporting Information (Table SA1).
We used the movement parameters from the iSSA to derive a selection-free movement kernel for turtles in each behavioural state.
As expected, turtles in the ARS state exhibited short step lengths and uniform turning angles, while turtles in the transit state exhibited longer step lengths and turning angles more concentrated around 0 ( Figure 4). ARS steps (in meters) were gamma distributed with shape = 1.05 and scale = 2728.42, and transiting steps (in meters) were distributed gamma with shape = 2.94 and scale = 10,080.28 (Figure 4a). Following from these distributions, expected step length during ARS was 2.9 km (SD = 2.8 km), and expected step length during transit was 29.6 km (SD = 17.3 km). ARS turning angles were von Mises distributed with mean = 0 (fixed) and concentration = 0.01, and transiting turning angles were von Mises distributed with mean = 0 (fixed) and concentration = 3.72 (Figure 4b).

DISCUSSION
While MPAs have consistently been shown to conserve ecosystem function and help promote sustainable fisheries (Roberts et al., 2018), gauging their effectiveness for migratory species is challenging. MPA networks that are designed to incorporate essential habitats used by TA B L E 1 Summary of parameter estimates included in the model and number of total locations per model incorporated in the analysis. Boxes with a (+) highlighted in dark grey indicate positive selection and therefore preference relative to unprotected waters. Boxes with a (−) highlighted in light grey indicate avoidance relative to unprotected waters. Boxes with a 0 indicate no selection/preference F I G U R E 2 Estimates of log-RSS for MPAs. RSS indicates how much more likely a step is to end in protected waters (multi-use or no-take) versus unprotected waters. The natural logarithm of RSS is the scale estimated by the iSSA, with a log-RSS of 0 (dashed black line) indicating no preference versus unprotected zones. Bars show 95% CIs for the log-RSS estimate. Turtles showed positive selection for multi-use zones in both behavioural states, whereas they showed no selection for no-take zones during either state (confidence intervals overlap 0) mobile species, such as spawning, nursing or foraging grounds, could fill a crucial role in maintaining population persistence. Therefore, MPAs targeting such protection should be systematically evaluated to quantify how migratory species interact with these boundaries.
Step selection analysis is a flexible statistical tool that can expand our understanding of animal movement and behaviour at fine spatial and Pacific (Kobayashi et al., 2008), suggesting these animals often seek out productive habitats while foraging. One possible explanation for the reverse outcome from our model could be due to the high degree of uncertainty associated with interpreting chlorophyll as a model parameter across large spatial or temporal scales (Brewin et al., 2017).
Monthly chlorophyll composites can contain a high percentage of cloud cover, large gaps or inflated values in nearshore environments due to bottom contamination (Blakey et al., 2016). Our study area is predominantly shallow (median depth ∼ 14 m) so it is likely that high chlorophyll signals are mixed with seagrass, algae or other substances in the water column. Furthermore, our data are biased towards loggerheads (68%), who forage primarily on invertebrates and could therefore avoid dense areas of seagrass or potentially high chlorophyll areas.
Designing MPAs will naturally involve compromise among competing objectives, specifically where no-take areas are concerned. If an objective is to partially or completely encompass the migratory corridors and/or foraging grounds of marine megafauna, precise placement of no-take areas is key. Our results did not show that turtles select for no-take areas, and there may be several explanations as to why. First, no-take areas in our study region are rare compared to multi-use zones, so our iSSA has a small sample size for both used and available locations falling within no-take zones, and the resulting parameter estimates thus have a large amount of uncertainty. Second, there may be an issue of location resolution. No-take areas in our study region are also relatively small, and Argos telemetry data are known to have low spatial precision. While the SSM largely accounts for this, there is still uncertainty in each location and the resulting steps, potentially making F I G U R E 3 Estimates of log-RSS for depth (a) and chlorophyll (b) for turtles in each behavioural state (ARS) = blue, transit = orange). RSS indicates how many times more likely a step is to end in the given depth (along the x-axis) versus the mean for that covariate (mean depth = -18 m; mean chlorophyll = 1.2 mg/m 3 ). The natural logarithm of RSS is the scale estimated by the iSSA, with a log-RSS of 0 (dashed black line) indicating no preference versus the mean. Coloured dashed lines show 95% confidence envelope for log-RSS estimates. Turtles showed positive selection for shallower depths during ARS, but no selection for depth during transit (a). Turtles showed avoidance of the highest chlorophyll concentrations during ARS, but no selection for chlorophyll during transit (b) F I G U R E 4 Selection-free movement kernels estimated by the iSSA by behavioural state. During ARS, turtles exhibited relatively short step lengths (a) and uniform turning angles (b). During transiting movements, turtles exhibited relatively long step lengths (a) and very directed movements (b; i.e. relative turns concentrated around 0) it difficult to detect selection for small sites (Fortin et al., 2005). Third, sea turtles may truly not select no-take areas, suggesting these zones may be poorly designed for these species. For example, sea turtles may not be able to recognize the boundaries of these protected areas since there can still be boat traffic, substantial human presence, bait fishing and, occasionally, catch-and-release trolling (U.S. Department of Commerce, 1996). Nevertheless, results of this analysis could be integrated into future rezoning considerations for sea turtles and other marine megafauna.
Our analysis does not allow us to address how or why turtles select for protected areas. The finding that sea turtles select multi-use areas during ARS and transiting implies that they can somehow perceive them, but this could be interpreted several ways. The protections in these areas might improve the overall health of the ecosystem and quality of the habitat within them, and turtles are simply selecting the highest quality habitat. Alternatively, turtles might be avoiding the disturbance outside of the MPAs, so our findings could be interpreted instead as an avoidance of unprotected habitat. Further investigation into these behaviours at a finer spatial and temporal scales -for example using passive acoustic receivers and/or GPS tags -would help to shed light on some of these important questions.
Previous efforts to inform resource management with animal tracking data have ensured that the developed models can account for significant variation at numerous scales (temporal, individual, etc.) and are accessible to managers, stakeholders and policy makers (Hays et al., 2019). While planning for effective MPAs requires the consideration of many unique components, satellite tracking studies have often had a significant influence on MPA design and placement (Lea et al., 2016;.

CONCLUSION
To design the most effective MPAs and conservation strategies for threatened and endangered marine species, a comprehensive understanding of animal movement patterns and habitat selection is necessary. Previous research has identified various techniques, including hotspot analyses and kernel density estimation, that could reveal areas of 'high use' within a migratory species' home range; this information can be used to better focus management efforts to target protections for the most vulnerable species (Lascelles, Langham, Ronconi, & Reid, 2012). However, the use of these more traditional methods often only provide static summaries of habitat-use patterns or limit/exclude environmental covariates that are critical to discerning species distribution and behaviour. It is now possible, although more complex, to incorporate behavioural data into spatial models for conservation planning purposes (Ashe, Noren, & Williams, 2010). Our spatially explicit model is an alternative method to evaluating MPA boundaries and zoning arrangements for migratory species. Our approach eliminates tagging site bias and provides behavioural context with movement patterns. This method may be adapted to any species with satellitetracking data and inclusive of a wider range of covariates. Our results emphasize the improvement in model predictions when satellitetracking data are integrated into the design and management of MPAs.