Telemetry‐validated nitrogen stable isotope clocks identify ocean‐to‐estuarine habitat shifts in mobile organisms

Throughout their life history, many animals transition among heterogeneous environments to facilitate behaviours such as reproduction, foraging and predator avoidance. The dynamic environmental and biological conditions experienced by mobile species are integrated in the chemical composition of their tissues, providing retrospective insight into movement. Here, we present a unique application of nitrogen stable isotope clocks (‘isotopic clocks’), which integrate tissue turnover rates, consumer stable isotope ratios and habitat‐specific isotope baselines to predict time‐since‐immigration and the timing of habitat shifts in a migratory species. Nitrogen stable isotope values of blood plasma collected from juvenile sand tiger sharks Carcharias taurus, a species known to undertake seasonal movements between ocean and estuarine environments, were used to derive estimates of time‐since‐immigration and the timing of seasonal habitat shifts undertaken by this species. Nitrogen isotopic clocks estimated for 65 juvenile sand tiger sharks sampled across 6 years indicated that individual sharks predominantly arrived to estuarine habitats between June and July, with some individuals arriving as early as mid‐May. These estimates were validated by comparing isotope‐derived estuarine arrival times with those from acoustically tracked individuals. The median estuarine arrival day estimates from our isotopic approach aligned with estimates from acoustic telemetry for each sampling population. Sensitivity analyses indicated that isotopically inferred time‐since‐immigration and estuarine arrival estimates were robust to variation in isotopic turnover rate and diet tissue discrimination factors under multiple modelling scenarios. This suggests that parameterization of the nitrogen isotopic clock provides reliable estimates of time‐since‐immigration and day of arrival into new habitats if isotopic variation exists between origin and new locations. Our study presents a unique application of telemetry‐validated isotope clocks to derive retrospective estimates of time‐since‐immigration and timing of habitat shifts for animals that seasonally traverse heterogeneous environments. This approach can be readily applied across many temporal and spatial scales, and to other species and ecosystems, to facilitate rapid assessment of changes in animal habitat use and broader ecosystem structure.


| INTRODUC TI ON
The physical, chemical and biological properties of environments through which animals move are often integrated in the biochemical composition of their tissues, either through diet (i.e. catabolized to organic macromolecules) or accretion into hard structures (e.g. bound within otoliths, teeth and bone). The potential of biochemical markers to record migration dynamics is great (e.g. Graham et al., 2010;Trueman et al., 2019), but application lags behind that of conventional tracking approaches (e.g. biotelemetry, Graham et al., 2010;Hobson et al., 2012;Madigan et al., 2014;Moore et al., 2016;Vander Zanden et al., 2018;Trueman et al., 2019). Naturally occurring stable isotope ratios form the methodological basis of most studies because regional variability in the isotope ratios of baseline organisms (defined as low trophic level prey species with high site fidelity), which is largely driven by local environmental (e.g. precipitation, temperature, salinity) and biogeochemical processes (e.g. nitrate vs. nitrite metabolism or C 3 vs. C 4 photosynthesis), is integrated into the tissues of higher organisms that readily migrate (Madigan et al., 2014, in press;Moore et al., 2016). In most cases, the tissues of an animal that has recently migrated from one region to another, assuming both are isotopically distinct, will reflect an admixture of origin and new prey isotope baselines. This is because the individual has not yet reached an isotope 'steady state' with the new system (Madigan et al., 2014), provided that sampling has occurred within the time window of isotope turnover for a given element (Thomas & Crowther, 2015;Vander Zanden et al., 2015). From examination of an individual's 'isotopic clock' (Klaassen et al., 2010), knowledge of these facets can be used to estimate time-since-immigration, and thus the retrospective timing of individual movement into a new system (Madigan et al., 2014(Madigan et al., , 2017Moore et al., 2016). Defining the timing of animal migration is important to establish temporal variability in food web dynamics, as well as to encourage cooperation between management boundaries over which individuals may move (Madigan et al., 2021).
Application of an isotopic clock has provided new insight into the migration dynamics of Pacific Bluefin Tuna Thunnus orientalis across the Pacific Ocean (Madigan et al., 2014). In Madigan et al. (2014), timesince-immigration and approximate arrival times were estimated for individuals moving from the isotopically distinct, nutrient-poor waters of Eastern Japan, with a lower baseline nitrogen isotope ratio (δ 15 N), into the highly productive, nitrate-dominated waters of the Eastern Pacific Ocean, with a higher baseline δ 15 N. The authors were able to define residents from migrants and evaluate whether migration patterns varied across ontogeny. Isotopic clocks have also been applied across marine-brackish-freshwater gradients; for example, Moore et al. (2016) combined carbon, sulphur and nitrogen (δ 13 C, δ 34 S and δ 15 N) isotopic clocks to predict time-since-immigration for juvenile life-history stages of five Pacific salmonids within the Skeena River watershed, British Columbia. This approach was used to confirm estuaries as critical stopover habitats for life-history stages that were too small to track using conventional approaches (Moore et al., 2016).
Despite the promise of isotopic clocks to discern important information on animal migrations and habitat use, there have been few applications to other species and study systems, especially at smaller spatial and temporal scales. Similarly, a critical comparison between isotopederived migration information with that of telemetry approaches (e.g. passive acoustic vs. isotopic clock) has yet to be undertaken.
We present results from both nitrogen isotopic clocks (δ 15 N) and acoustic telemetry to illustrate how isotopic clocks can be used to estimate time-since-immigration and seasonal movement between heterogeneous marine habitats. Juvenile sand tiger sharks Carcharias taurus are used as a case study species as this species performs seasonal ocean-to-estuarine habitat shifts (Kneebone et al., 2012(Kneebone et al., , 2014. Our application to juvenile sharks illustrates that stable isotope clocks are applicable at much finer temporal and spatial scales than documented in previous studies (Madigan et al., 2014;Moore et al., 2016).
Our study provides a clear methodological framework that can be applied to other animals that traverse chemically distinct interfaces, to estimate time-since-immigration and infer the timing of habitat shifts, collectively informing the broader ecology of species and ecosystems.

| Study location and sampling of juvenile sand tiger sharks
The Great South Bay (GSB), New York, US is a highly productive temperate lagoon that provides critical stopover habitat and food resources for many juvenile life-history stages of marine finfish, especially during the summer and early fall (Olin et al., 2019;Tinoco, 2017).
Anecdotal evidence dating back to the early 20th century proposed that GSB has long served as a potential summer nursery ground for juvenile sand tiger sharks Carcharias taurus (Thorne, 1928). Juvenile sand tiger sharks are known to migrate considerable distances between southern wintering grounds, extending from Florida to Northern Carolina, and summer nursery grounds in more northern, temperate embayments (Kneebone et al., 2014; Figure 1). In coastal New York, approach can be readily applied across many temporal and spatial scales, and to other species and ecosystems, to facilitate rapid assessment of changes in animal habitat use and broader ecosystem structure.

K E Y W O R D S
Carcharias taurus, ecosystem connectivity, habitat use, heterogeneous habitats, movement ecology, shark their migration dynamics remain speculative, but offer an opportunity to develop an isotopic clock model to discern the dynamics of ocean-estuarine movements. Juvenile sand tiger sharks were captured between May and September in consecutive years from 2014 to 2019 in GSB (Figure 1) using rod and reel angling and drum lines. Some individuals were recaptured multiple times throughout the six sampling years, and these data were treated as independent fish for the purposes of the isotopic clock analysis. A full description of animal sampling is given in the Supporting Information (Appendix S1).

| Isotopic clock and statistical resampling procedure
Isotopic clocks are defined by an equation from Klaassen et al. (2010) that describes time as a function of the difference in ecosystem baseline isotope ratios and the rate of turnover for that isotope in an organism: where t = time (in days) an individual has been resident in the new system, δ 0 = the assumed isotopic composition of an animal at isotopic steady state with their origin location, δ f is the assumed isotopic composition of the animal at isotopic steady state with the new location, δ t is the isotopic composition of an individual at the time of sampling and λ is the tissue-specific isotopic turnover rate ( Figure 2). This equation can be solved by entering a value for each of the input parameters to provide a single estimate. Alternatively, t can be viewed as a distribution of estimates drawn from distributions of input parameters that represents a more realistic characterization of uncertainty (see Heady & Moore, 2013;Klaassen et al., 2010). Thus, t i is 1:10,000 draws from normal distributions of δ 0i , δ fi and λ i based on the mean (µ) and standard deviation (σ; Figure 2). The first 1,000 successful estimates of t i are then used to generate distributions of time-since-immigration, t i (in days) and day of arrival with the following equation: To estimate population-level estimates of t i and arrival day, the distributions for all individuals can be combined for each year (e.g. Madigan et al., 2014;Madigan et al., 2021). We further explore uncertainty in t i estimates by conducting sensitivity analyses related to tissue turnover rates and diet tissue discrimination factor/trophic enrichment factor (herein DTDF) as described below. Diet tissue (2) Arrival Day = Sampling Day − t i .

F I G U R E 1
Proposed migration route of juvenile sand tiger sharks into Great South Bay, New York from the New York Bight. Top inset depicts the isotopic biplot for juvenile sand tiger sharks captured in GSB by year in relation to δ 0 and δ f parameters used in the isotopic clock model. Bottom inset represents the proposed seasonal wintering ground of juvenile sand tiger sharks; δ 0 were parameterized based on the isotopic composition of blood plasma sampled from juvenile blacktip sharks from the First Coast, Jacksonville, FL. Source: The Great South Bay and NE US Seaboard, Google Earth, Accessed 9 January 2020 discrimination factors are defined as the stepwise isotopic enrichment/depletion in a consumer relative to its prey (Caut et al., 2009;Hussey et al., 2010). First, we describe the development of parameter estimates that populate-Equation 1.

| Model parameterization for juvenile sand tiger sharks
We used isotopic clocks to generate time-since-immigration estimates for juvenile sand tiger sharks captured in GSB. Models were parameterized using a single elemental isotope ratio, δ 15 N, measured in shark blood plasma, for two reasons. First, the high terrestrial nutrient input into estuaries often results in elevated 15 N compared to primary producers of continental shelf ecosystems. For GSB specifically, high nitrogen loading from local agriculture (e.g. duck farms, Gobler et al., 2019) provides sufficient isotopic distinctiveness relative to the coastal ocean, which is required to effectively parameterize isotopic clock models (Madigan et al., 2014, in press).
In coastal marine settings, a similar range in δ 13 C between baselines is certainly possible, but high variability of planktonic and benthic primary producers that utilize both C 3 and C 4 photosynthetic pathways may be problematic for parameterizing accurate isotopic baselines. Second, the use of δ 13 C isotopic clocks in this study system may be further compounded by uncertainty surrounding the lipid dynamics of shark blood plasma . Given the potential for lipid effects, δ 13 C was assumed to be less reliable compared to δ 15 N for the parameterization of isotopic clocks. Isotopic clocks were implemented for animals captured from 2014 to 2019 ( Figure 1; Table 1). The isotopic turnover rate used for 15 N of shark blood plasma (expressed in day −1 ), λ, was 0.0165 ± 0.003 day −1 based on estimates derived from leopard sharks Triakis semifasciata F I G U R E 2 Flow diagram illustrating parameterization and implementation of an isotopic clock approach to calculate time-sinceimmigration and arrival day for highly migratory animals. Note that δ 0 and δ f parameters can be derived from prey-based sampling + DTDFs or via sampling functionally equivalent species assumed to be at isotopic steady state with the origin or new locations held in experimental conditions (Kim, del Rio, et al., 2012). We used the isotopic turnover rate derived from this study as it represents one of the only published estimates for shark blood plasma that is coupled with estimates of diet tissue discrimination factors that can be used to parameterize δ f .
Isotopic baselines or functionally equivalent species can be used to define the hypothetical isotopic composition of juvenile sand tiger sharks at relative isotopic 'steady state' with the origin and new ecosystem (Madigan et al., 2021). These 'baselines' reflect the isotopic ratio of species that are assumed to inhabit the same trophic position, thus exhibit similar diet tissue discrimination from dominant source/s of primary production (Moore et al., 2016) or are based on potential prey items corrected for the DTDF (Heady & Moore, 2013;Madigan et al., 2014). In this study, we used existing isotope data from two independent datasets (Table 1; Appendix S1 Table S1) to derive δ 0 and δ f .
Existing biotelemetry data suggest that the coastal waters of North Carolina and Florida serve as seasonal wintering grounds for juvenile sand tiger sharks (Kneebone et al., 2014). Thus, we assumed δ 0 (i.e. reflecting the animal at isotopic steady state with the origin location) would reflect juvenile sand tiger sharks at isotopic steady state with isotopic prey baselines reflecting southern wintering grounds.
Due to a lack of isotopic data from juvenile sand tiger sharks in this region, we parameterized (δ 0 ) based on literature-derived δ 15 N of blood plasma measured for a functionally equivalent species known to exhibit similar dietary habitats and migration dynamics (based on published diet data, Cortés, 1999), the blacktip shark Carcharhinis limbatus (δ 15 N = 12.7‰, Morgan et al., 2020). Owing to the potential relationship between body size and δ 15 N as commonly observed in marine fishes (Jennings et al., 2002;Galvan et al., 2010), literature-derived values were based on small-bodied individuals that exhibited a similar size range (fork length = 68-155 cm) to juvenile sand tiger sharks sampled in this study (fork length = 89-158 cm). Blacktip sharks were sampled from a central portion of the assumed wintering grounds, the First Coast, Jacksonville, FL (Morgan et al., 2020). We assigned a conservative standard deviation of ±1.0‰, which is higher than that reported by Morgan et al. (2020, SD = 0.83‰) to ensure our resampling procedure adequately accounted for potential variation in isotopic baselines.
Estuarine baselines were derived from the isotopic composition of white muscle tissue of small-bodied secondary and tertiary consumers (Appendix S1: Table S1), including weakfish Cynoscion regalis (n = 19, TL = 6.5 cm-26.5 cm), summer flounder Paralichthys dentatus (n = 31, TL = 11.4 cm-38.5 cm) and winter flounder Pseudopleuronectes americanus (n = 31, TL = 7.2 cm-21.0 cm). These species were selected because they likely serve as potential prey items for juvenile sand tigers in estuarine food webs, given the latter's broad dietary preferences noted elsewhere (Cortés, 1999). Even if juvenile sand tiger sharks do not directly prey on these species, their isotopic composition represents that of low to intermediate consumers within the estuarine food web and is assumed synonymous to juvenile sand tiger shark diet. Additionally, our sensitivity analyses (see below) ensured that an extremely broad isotopic range of potential prey could be accounted for in our modelling framework. The mean isotopic composition of all individuals combined (n = 81, δ 15 N = 13.4‰) was corrected based on published DTDFs for blood plasma (Kim et al., 2012, see  This represents one of the only isotopic turnover rates published for blood plasma that is derived from controlled feeding experiment. This study is also the source of DTDFs used for the parameterization of δ f Kim, del Rio, et al. (2012) Δ 15 N = 3.2‰ because this value represents the mid-point of published estimates for leopard sharks where individuals were fed single diets of squid (DTDF = 2.2‰) and tilapia (DTDF = 4.2‰) in a laboratory setting . We assume that juvenile sand tiger sharks in natural systems feed across multiple prey types and are therefore more likely to exhibit an intermediate DTDF reflecting a mixed diet. We assigned a conservative standard deviation of ±1.0‰ to δ f , which ensured that variation associated with DTDFs and isotopic baselines was accounted for (Phillips et al., 2014). Furthermore, the effects of uncertainty on estimated time-since-immigration and estuarine arrival were tested across several sensitivity models (see below).

| Sensitivity analysis
Isotopic turnover rates (λ) and diet tissue discrimination factors It is therefore recommended that such models are examined for sensitivity to the use of variable DTDFs (Olin et al., 2013;Phillips et al., 2014).

| Acoustic telemetry
We examined the accuracy of our primary isotopic clock model by comparing model predictions to validated passive acoustic telemetry data, for each sampling year. Specifically, time-of-first detection was calculated from an existing passive acoustic telemetry array, which comprised between 19 and 28 VR2W (69 KHz) acoustic hydrophones (Vemco Inc). These were deployed from April to November between 2014 and 2019 (Appendix S1: Figure S1, Table S2).

| RE SULTS
Blood plasma/serum was collected from 66 juvenile sand tiger sharks between 2014 and 2019 (Table 2); this included recaptures of 14 individuals across multiple sampling years. May to September, Appendix S1: Table S3), overlapped considerably with probable arrival days inferred from the nitrogen isotope clocks.

| Primary model
In all cases, the range of entry times captured the median (i.e. 50%) arrival day for the sampled population within each year (Figure 4).

| Sensitivity analysis
The constant turnover rate isotopic clock model provided median time-since-immigration and estuarine arrival days for population-level estimates that were within ≤2 days of the primary model where isotopic turnover was assigned a SD of 0.003 (Figure 4; Table 3; Appendix S1: Table S4). There were small differences in the range of possible time-since-immigration and arrival day estimates between the primary and constant turnover rate model as inferred through percentiles (Figure 4; Table 3; Appendix S1: Figure S2, Table S4). Across all sampling years, assuming a constant isotopic turnover rate generally reduced the upper 97.5th percentile by ≤18 days; there was a negligible effect on the lower 2.5th percentile ( Figure 5, Table 3). Time-since-immigration (t i ) and estuarine arrival estimates were somewhat sensitive to the use of variable DTDFs assigned to δ f , especially when the assigned DTDF was lower than the primary model (Table S5). Assigning a lower DTDF of 2.2‰ (δ f = 15.6 ± 1.0‰) resulted in longer time-since-immigration estimates and earlier estuarine arrival estimates by ≤26 days relative to the primary model, which is relatively F I G U R E 4 Density plots highlighting the sensitivity of time-since-immigration (t i ) and estuarine arrival day across a series of model scenarios, relative to the primary model (green): a constant turnover rate model (pink), and two DTDF sensitivity models where δ f = 15.4 ± 1.0‰ (red) and 17.4 ± 1.0‰ (blue), respectively. Models were run to provide population-level estimates of time-since-immigration and estuarine arrival for each sampling year  small compared to the total time-since-immigration time (Figure 3;

F I G U R E 3 (a)
Appendix S1: Table S5). This model also increased the uncertainty surrounding the median time-since-immigration and arrival times for each year as inferred through percentiles (Table 3; Appendix S1: Figures S2 and S3; Table S5). Assignment of a higher DTDF of 4.2‰ ± 1.0 (δ f = 17.6 ± 1.0‰) resulted in more agreeable median time-sinceimmigration and estuarine arrival estimates (≤15 days) relative to the primary model ( Figure 4; Table 3; Appendix S1: Figures S2 and S3, Table S6). Model uncertainty was also lower relative to the primary, constant turnover and lower bound DTDF models (Figure 4; Table 3; Appendix S1: Figure S2 and S3, Tables S4, S5 and S6).

| D ISCUSS I ON
Nitrogen stable isotope clocks discerned the timing of ocean-estuarine habitat shifts and provided time-since-immigration estimates for a highly migratory marine species. We demonstrated that the timing of habitat shifts in juvenile sand tiger sharks, as inferred through stable nitrogen isotopes (i.e. estuarine arrival days), was consistent with passive acoustic telemetry, suggesting that stable isotope clocks can provide complementary information to traditional tracking methods in near-shore settings. This approach has application to a broad range of organisms that seasonally traverse heterogeneous environments in terrestrial, freshwater and marine systems (e.g. crustaceans, fishes, turtles and mammals). Below we discuss the various model considerations and further validation that is required to promote the increased and broad-scale use of this approach in future migratory studies.

| Model considerations and validation
The ability to rapidly calculate the timing of habitat shifts at relatively low cost has the potential to greatly expand the study of animal migration, and offers an attractive means for studying species that may be too small to tag using traditional approaches (Moore et al., 2016).
Here, we provide the first comparison between stable isotope-derived arrival days with those inferred through passive acoustic telemetry.

F I G U R E 5
Density distributions highlighting annual entry times for juvenile sand tiger sharks modelled using isotopic clocks (density plots) and range of first detection times gathered from individuals detected on a long-term passive acoustic telemetry array. Dashed lines represent the full range of acoustic-derived entry times and solid lines represent the median estimates for each year. The bottom right panel illustrates scientists from the New York Aquarium conducting routine maintenance on their passive acoustic telemetry array. Photo: Julie Larsen Maher © WCS We observed that the most probable isotope-derived arrival days (at the population level) fell within the range of first detections estimated through acoustic telemetry. This provides strong support that isotopic clock models can accurately predict the arrival periods of juvenile sand tiger sharks into the summer nursery grounds of GSB.
A current limitation of stable isotope-based modelling is the ambiguity surrounding diet tissue discrimination factors and tissue turnover rates for commonly analysed tissues (Gannes et al., 1997;Shipley & Matich, 2020). In many cases, tissue turnover rates and DTDFs are unavailable for the study species, necessitating the use of literature-derived values from taxonomic relatives (Healy et al., 2018 estimates. This is because the differences in isotopic baseline composition between origin and new locations, in most cases, are greater than the isotopic variation associated with individual-level foraging (Madigan et al., 2021). Therefore, assigning appropriate isotopic baselines in isotopic studies -remains a -challenge owing to the dynamic nature of baseline composition. Future studies applying isotope clocks should consider these challenges and aim to incorporate robust evaluations of animal ecologies (i.e. extent of individual-level variation in diet) and isotopic baselines of varying ecosystems (and inherent variation in baselines), as well as the application of multiple isotopes or expansion of additional tracers into the modelling framework (see Heady & Moore, 2013).

| Application to species and ecosystems
Understanding the timing and duration of seasonal habitat use is increasingly needed as researchers seek to understand the effect of climate change on species phenology and the resulting response by ecosystems (Staudinger et al., 2013;Weiskopf et al., 2020). The current research was conducted in the Great South Bay, NY, an ecosystem that is driven by the timing of a suite of predator species that occupy the ecosystem seasonally and include sand tiger sharks, summer flounder, black seas bass Centropristis striata, and striped bass Morone saxatilis, among others (Nuttall et al., 2011;Olin et al., 2019).
Here, the ability to rapidly identify the timing of juvenile shark movements into this system is fundamental to local conservation initiatives; this is because estuarine systems encompass important nursery habitat for juvenile sharks during certain times of the year. Therefore, knowledge of their arrival and time-since-immigration is imperative for providing regional management of their populations and protection of key habitat types.
To fully appreciate the changes to predator-prey relationships and ecosystem dynamics, methods are needed to quantify the timing of movement between ocean and estuarine environments. The case study presented herein provides confirmation that isotopic clocks can be applied to estimate the timing and duration of habitat shifts in a seasonal predator in a coastal ecosystem. Using confirmation through acoustic telemetry, we illustrate the efficacy of isotope clocks to measure shifts between heterogeneous habitats at finer temporal and spatial scales than previous applications of this approach. Overall, this work emphasizes that a broader application across animal taxa and ecosystems is possible.

ACK N OWLED G EM ENTS
The authors would like to thank C.