Seasonal patterns in stable isotope and fatty acid profiles of southern stingrays (Hypanus americana) at Stingray City Sandbar, Grand Cayman

Ecotourism opportunities in the marine environment often rely heavily on provisioning to ensure the viewing of cryptic species by the public. However, intentional feeding of wildlife can impact numerous aspects of an animals’ behavior and ecology. Southern stingrays (Hypanus americana) provisioned at Stingray City Sandbar (SCS) in Grand Cayman have altered diel activity patterns and decreased measures of health. This study looked at seasonal changes in stable isotope (SI) and fatty acid (FA) profiles of provisioned stingrays at SCS. Plasma δ15N was higher in male stingrays (11.86 ± 1.71‰) compared to females (10.70 ± 1.71‰). Lower values for δ15N in males and females were measured in October during low tourist season, suggesting stingrays may be forced to rely on native prey items to supplement the decreased amount of provisioned squid available during this time. Plasma FA profiles were significantly different between sexes and across sampling time points, with FAs 22:6n3, 16:0, 20:5n3, 18:1n3C, 18:0 and 18:1n9T contributing to dissimilarity scores between groups. Dietary FAs primarily contributed to differences between males and females lending further evidence to differences in foraging patterns at SCS, likely due to intraspecific competition. Further, canonical analysis of principal coordinates (CAP) analysis of FA profiles suggest similar diets during peak tourist season and differences in diet between males and females during the low season. This study demonstrates alterations in feeding ecology in stingrays at SCS which is of critical importance for effective management of the SCS aggregation.


Results
The female stingrays randomly selected for the study averaged 107.2 cm disc width (DW), while the males were considerably smaller, as expected with this sexually dimorphic species (53.8 cm DW). All animals had a PIT tag at the initial capture in January (10 males, 10 females) and once it was discovered that re-capture rate in April was poor for males (20% re-capture rate of males [2 individuals], 18 females), animals were sampled opportunistically in the subsequent months of July (7 males, 18 females), and October (7 males, 17 females).
Males and females within sampling month were distinguished by FA profiles (CAP p = 0.001). Fatty acids identified by SIMPER routines were overlaid as vectors on the CAP ordination plot (Fig. 3).

Discussion
Few studies have explored stable isotopes in species belonging to the Family Dasyatidae. Nitrogen isotope compositions in the muscle of brown stingrays Dasyatis lata in Kane' ohe Bay, Hawai'i (range δ 15 N 10.4 to 11.9‰, δ 13 C −13.5 to −11.9‰ 45 ) were similar to those found in southern stingrays in our study (range δ 15 N 9.7 to 12.6‰, δ 13 C -18.2 to -16.7‰ in plasma). However, muscle sampled from southern stingrays in Belize had lower δ 15 N values (mean δ 15 N 7.4‰ and δ 13 C -8.5‰ 46 ), which the authors suggested may have resulted from a diet of softbodied invertebrates and bivalves rather than crabs and teleost fish typically consumed by southern stingrays in other locations 37 . δ 13 C values for southern stingrays at SCS were lower compared to values reported for brown Table 1. Mean (± SD) δ 15 N and δ 13 C (‰) measured in the plasma and red blood cells of male and female southern stingrays sampled in January, April, July, and October at Stingray City Sandbar. Sample sizes in parentheses.   www.nature.com/scientificreports/ stingrays in Hawaii and southern stingrays in Belize and differences likely reflect different habitats and sources of carbon at the base of each food web 16 . The heavy nitrogen isotope ( 15 N) increases in concentration in a stepwise fashion with increases in trophic level. Conversely, δ 13 C is conserved through trophic systems, and δ 13 C values are predominantly used to determine the source of carbon from primary producers 16,17,47 . Fractionation values differ for different species and tissues sampled and few studies have measured fractionation values in elasmobranchs. Only a single study has examined fractionation values in batoid rays. Controlled feeding trials conducted on smallnose fanskates (Sympterygia bonapartii) resulted in fractionation values of 1.7 ± 0.1‰ for δ 13 C and 2.5 ± 0.2‰ for δ 15 N in whole blood 48 . Plasma of male stingrays at SCS was found to be higher in δ 15 N compared to females (11.9 ± 1.7‰ vs. 10.7 ± 1.7‰) in our study and males had significantly higher δ 15 N than females year-round, except in April. Lack of a difference between the sexes in April likely reflects the small number of males that were captured (n = 2). Similar trends were observed for δ 15 N values in red blood cells, where males had higher δ 15 N compared to females in all months with the exception of April. Higher δ 15 N in males may be indicative of males feeding at a higher trophic level (although fractionation values do not indicate a full trophic level increase) than females, nutritional stress in males, differences in growth rates between males and females and/or differences in diet quality, both of which can influence δ 15 N discrimination values 49 . Squid offered at the SCS was not sourced locally and therefore may have δ 15 N values that do not reflect the regional food web. Based on stomach content analyses [36][37][38]50,51 , stable isotope values of potential prey of southern stingrays from the US, Bahamas, Belize, and Cuba, were extracted from the literature (Table S2) and plotted with squid values (Illex and Loligo sp.) from the North Pacific and North Atlantic (Fig. 4) to understand what may be driving these differences in nitrogen signatures between males and females at SCS.  www.nature.com/scientificreports/ Higher δ 15 N in the plasma and red cells of males in this study suggests males may be receiving more provisioned squid in their diet than females. This is unexpected given the differences in movement patterns of male and female stingrays at the SCS 6 and despite a clear preference for the provisioning of larger females by tour operators (L. Hoopes, personal observation). Corcoran et al. 6 noted that males moved actively along reef edges at night while females remained inactive in close proximity to SCS. They suggest that larger females, who were targeted by tour operators, were receiving greater quantities of squid and out-competing males for access to provisioned squid, thus requiring them to engage in nocturnal foraging 6 . Similar competition between sexes for tourist-derived food was observed in Australia where large female Dasyatis stingrays behaved aggressively towards smaller males 52 . If males were foraging nocturnally to supplement their diet with natural prey items, their diet would likely be comprised of crustaceans (portunid crabs, penaeid shrimp) and teleost fishes (Families Labridae, Gobiidae, Scaridae) which accounted for 87% of stomach contents of southern stingrays in the Bahamas 37 . However, δ 15 N values of native prey from the literature range broadly from 4 to 8‰, at least one trophic level below the values sampled in males in this study (Fig. 4). Potential native prey was not sampled during this study and SI values for invertebrates and teleosts within the Cayman Islands could not be found in the literature. While prey values from other locations may not accurately reflect baseline δ 15 N within SCS and surrounding reef systems, male stingrays would need to be foraging on a native diet at an unrealistically high trophic level to reflect the blood nitrogen isotope composition measured in this study. Non-provisioned stingrays were not able to be sampled as controls in this study and we cannot eliminate the possibility that differences in δ 15 N between males and females were an artifact of variables not related to provisioned feeding.
Elevated δ 15 N signatures can be indicative of nutritional stress in the form of caloric restriction, fasting, or starvation [53][54][55][56] . During periods of nutritional stress, animals may catabolize their own 15 N-rich tissues (essentially feeding at a higher trophic level) causing an increased excretion of 15 N-depleted nitrogeneous wastes, resulting in tissues elevated in 15 N 17 , including short-term increases in blood δ 15 N 53,57,58 . Studies attempting to measure higher δ 15 N during nutritional stress have been equivocal (e.g., 55,59 ). Logan and Lutcavage 60 measured temporal increases in blood δ 15 N in skates (Leucoraja sp.) who lost weight due to poor eating in the control group of a 20-d experimental turnover trial. Conversely, Wyatt et al. 61 noted a reduction in δ 15 N in the plasma of whale sharks (Rhincodon typus) fasting 151 days and suggested that urea may have an influence on tissue isotopes in calorically www.nature.com/scientificreports/ limited animals, although this individual later died and may have had other complicating health factors. Visual assessment of body condition of male stingrays during handling did not indicate starvation (L. Hoopes, personal observation) and higher δ 15 N in males due to severe caloric restriction or starvation seems unlikely. Differences in growth rates can affect stable isotope turnover and tissue discrimination and higher growth rates reduce nitrogen discrimination 62 . While male southern stingrays do exhibit faster growth rates than females in this sexually dimorphic species 39 , animals sampled in this study were within reported ranges for maximum disc width for adult males and females. This suggests that these animals had reached maturity were not likely actively growing and δ 15 N values would not be greatly affected by this variable.
Variation in δ 15 N may reflect differences in diet quality. In this context, diet quality is defined as the amino acid composition of a food source relative to the needs of a consumer 63 , but should also consider proximate composition (protein, fat, and carbohydrate balance) and the vitamins and trace minerals that are available across a range of prey items. Generally, carnivores feeding across a diversity of prey items would be expected to have high quality diets since they tend to consume prey with similar protein and amino acid compositions to their own tissues 61 . Higher quality diets result in dietary amino acids being routed directly to tissues with little processing and thus low isotopic discrimination. Conversely, lower quality diets (low protein with highly imbalanced amino acid compositions compared with consumer requirements) can inflate fractionation values of 15 N compared to high quality diets since consumers would need to synthesize scarce amino acids de novo from surplus amino acids 49,63 . As secondary consumers, squid (Ilex and Loligo sp.) have crude protein levels within the range of other teleosts (L. Hoopes, unpublished data), however, amino acid compositions are variable with notable differences in the proportion of glycine, proline, and arginine in squid compared to other marine prey 64 . It is unclear if the amino acid composition of a predominately squid diet meets stingray dietary requirements compared to a more varied diet including potential native prey at lower trophic positions, but it is unlikely that isotopic differences between males and females reflect greater de novo synthesis of certain amino acids. Regardless, the reduction in isotopic discrimination in higher order carnivores remains and area where further study is warranted 65 and future studies using compound specific amino acid stable isotope analysis may be useful in understanding whether diet quality differences exist between males and females feeding at SCS.
To evaluate the temporal scale in which isotopic values reflect diet, the turnover rates of those tissues must be known 61 . Isotopic turnover rate is the amount of time it takes for a given consumer tissue to incorporate the isotopic composition of the diet. Turnover rates vary depending on species, type of tissue being sampled, metabolic rate and animal size, and environmental variables like water temperature 61,66-70 . Assessing turnover time requires controlled feeding studies which are challenging to conduct, as they require long study periods and isotopically distinct diets 71 . Results of the few studies on turnover rates in elasmobranch species indicate that for plasma, the estimated days to 50% incorporation (half-life) of the new isotopic signature is approximately 35-60 days, although, whale sharks appear to have longer blood turnover rates (90-200 d) 48,60,61,72,73 . δ 15 N was significantly lower in October compared to other months in both red blood cells and plasma, indicating a probable shift in the diet of both males and females. October coincides with low numbers of tourists visiting the SCS (Port Authority of the Cayman Islands, www.cayma nport .com), presumably leading SCS stingrays to supplement their diets during this time feeding on a greater proportion of wild prey. Diet studies are lacking for southern stingrays from Grand Cayman, therefore using stable isotope values to infer consumption of specific prey must be interpreted cautiously. Lower δ 15 N signatures in invertebrate and teleost prey from the literature, compared to δ 15 N values in squid (Fig. 4) could suggest a greater reliance on native prey items during the low tourist season as SCS.
There were no differences in δ 13 C between males and females across tissues for any sampling month, with the exception of the month of July when plasma in males was slightly higher in δ 13 C (mean 0.12‰) than females. Carbon isotope compositions can become higher as a predator moves from feeding in an offshore food web to an inshore source 74 , however, these animals have shown high levels of site fidelity around the SCS site 6 and these detected differences are probably not biologically significant.
Fatty acid profiles for stingray plasma were significantly different across the year and between sexes and SIMPER results consistently identified dietary FAs (22:6n3, 18:1n9C, 20:5n3) contributing to these patterns. Fatty acid profiles have previously been examined in this population of stingrays 43 82,83 ), in various species of marine annelids (2.6 to 5.2% 84 ) and in conch species (2 to 5%, L. Hoopes unpublished data, 84 ).
While higher δ 15 N in males suggested greater consumption of provisioned squid compared to females, finer scale dietary changes were detected in FA profiles suggesting some scope for dietary differences (11-24%, SIM-PER) between males and females depending on season. Fatty acid profiles for males and females in October were different from stingrays sampled at other times of the year (Fig. 4) www.nature.com/scientificreports/ copepod-consuming species like crustaceans, bathypelagic squid, and fish in the diets of deepwater elasmobranch species 85,86 . While 18:1n7 can be metabolized from 16:1n7, Beckmann et al. 87,88 found this FA to be a reliable marker during controlled feeding trials with Port Jackson sharks (Heterodontus portusjacksoni). Monounsaturated FA 20:1n11 is typically found in high levels in cold water fish species 18,31 , and in lower levels in tropical and subtropical fish 89 , however, median levels of this FA were not significantly different across sampling time points. High levels of FA 20:4n6 (arachidonic acid, 13%) have been associated with Port Jackson sharks feeding on macro-algae consuming invertebrates (sea urchins and sea snails), while piked dogfish (Squalus megalops) known to feed exclusively on cephalopods had lower (4.6%) tissue levels of 20:4n6 90 . Conversely, Semeniuk et al. 42 noted higher levels of 20:4n6 in the plasma of native feeding southern stingrays (11.5%) compared to provisioned squid eating animals at SCS (5.2%). The incorporation of FAs from prey into the tissues of the consumers can be a complex process to understand in the elasmobranchs, however, the individual FAs driving these differences may suggest more reliance on native prey items in October compared to other months. October is considered the low season for tourists at SCS (Port Authority of the Cayman Islands, www.caymn aport .com) and if touristprovided squid are not as plentiful, foraging for wild prey may be increased. Conversely, FAs 20:5n3 and 16:0 were drivers of FA profiles in stingrays sampled in January and July (Fig. 3). Fatty acid 20:5n3 has been previously used as a biomarker of diatom-based food webs 33 and controlled feeding trials in Port Jackson sharks show it to be a more accurate indicator of diet in blood, compared to dietary FA 22:6n3 88 . Semeniuk et al. 42 noted high proportions of 20:5n3 in cold water elasmobranchs and stingrays at SCS fed a non-native diet of squid (25.5%). Fatty acids 20:5n3, 22:6n3, 18:1n9 and 18:2n6 have previously been identified as useful indicators of diet in the liver and muscle of deepwater chrondrichthyans 85 and in the liver, muscle, and blood of Greenland sharks (Somniosus microcephalus 33 ). January is considered the height of tourist season at SCS and visitor numbers to the sandbar reach yearly highs (Port Authority of the Cayman Islands, www.cayma nport .com). This may represent a time in which provisioning is at its peak and ostensibly, squid are plentiful for both male and female stingrays. Curiously, during the month of January, male and female FA profiles differed by 24% (SIMPER) suggesting there may be some finer scale seasonal dynamics at play within this population.
Although FA assimilation in sharks and rays appears to be similar to other taxa 87,88 , elasmobranchs rely on ketone bodies instead of lipids for energy metabolism and rely on lipid-rich livers for long-term energy storage 91 . The influence of caloric restriction and/or fasting on plasma fatty acids signatures in sharks is unknown. Wood et al. 92 noted decreases in non-esterified FAs in the plasma of dogfish sharks (Squalus acanthias) after feeding compared to fasted control animals, but they did not examine individual FAs and further research is warranted. The functions and fates of individual fatty acids during elasmobranch feeding and fasting may be taxa-specific and controlled feeding studies would further enhance the quantitative capacity of FA profiling in sharks and rays 93 .
Stable isotope data suggests a consistent, long-term pattern of provisioning was occurring at the SCS tourist site, however, FA profiles of SCS stingrays suggests small dietary differences between males and females during certain parts of the year, although natural temporal or spatial variations cannot be ruled out. Precise driving factors of these dietary differences remain unclear, but appear to be related seasonal differences in the quantity of squid available at the SCS and the nature of feeding on natural prey away from the SCS. Similar FA profiles between sexes in January suggest that both males and females are relying on provisioned squid during the height of tourist season, but when the number of tourists visiting SCS declines (October), both sexes may consume a larger proportion of natural prey items. Larger sample sizes for males as well as a SI and FA database of potential wild prey items would provide additional information for explaining feeding history and for predator-prey comparisons of profiles using Bayesian mixing models and quantitative fatty acid signature analysis 18 . Decades of provisioning have altered the feeding ecology of stingrays at SCS and considering the great economic value of the aggregation of stingrays, understanding the long-term effects of provisioning on the biology and overall health of the stingrays is of critical importance for well-informed and effective management of the SCS aggregation to ensure continued success of interaction activities.

Methods
This study was conducted with approval from the Cayman Islands Government Department of Environment. Animal handling and sampling was approved through the Research Committee at Georgia Aquarium and all handling procedures were conducted using guidelines established by the American Fisheries Society and American Society of Ichthyology and Herpetology. Study site. The SCS is a small, shallow sandbar (7800 m 2 ) with minimum depth of 0.5 m and sandy bottom located in the North Sound of Grand Cayman (19°22′34″N, 81°18′18″W) (Fig. S1). North of the sandbar is an adjacent patchwork of deeper reefs and habitat to the south is primarily seagrass beds. There are over one million visitors to the SCS site annually, where they interact with approximately 100 stingrays, 80% of which are large mature females (personal communication Cayman Islands Department of Environment, 39 ). Stingrays at SCS show very high fidelity to the feeding site, with females resident for an average of 10 years, and males resident for an average of 3 years 39 , demonstrating chronic provisioning of the animals.
Sampling and handling. Ten female and ten male stingrays ranging in disc width from 55 to 121 cm were chosen at random in January 2014 for this study. All animals had an existing passive integrated transponder (PIT) tag suggesting they were a known resident of the SCS population. Subsequent 4-day trips in April, July, and October of 2014 attempted to locate and re-sample these same 20 individuals, with mixed success. Between January and April, re-capture success of the original 10 females was 80%, while only 2 of the original 10 males were able to be located and re-sampled. Based on this low re-capture rate of males (20% re-capture in April), we opportunistically sampled males, and additional females, the remainder of the study period.

Scientific Reports
| (2020) 10:19753 | https://doi.org/10.1038/s41598-020-76858-w www.nature.com/scientificreports/ Stingrays were captured with hand-held dip nets and transferred into a seawater-filled pool inside a boat for examination and sample collection as part of a larger study. Once transferred onboard, handling time averaged 11 min and consisted of scanning for the presence of PIT tags, standard morphometrics, blood sampling, ultrasound assessment, and PIT tagging of any new individuals. Upon completion, stingrays were placed back in the dip net and returned to the water. Visual observation of released stingrays indicated that they typically resumed feeding within several minutes of being returned to the water from nearby tourist feeding interactions.
Blood was collected either ventrally from the caudal tail vein (21 g × 1.5 in. needles) or dorsally from a pectoral fin vessel (21 g or 23 g × 1 in. butterfly needle) and then transferred into lithium heparin vacutainers. No difference in blood isotope or fatty acid profiles were expected between venipuncture sites based on similarities in blood chemistry values from paired samples (L. Hoopes, unpublished data). Blood volume ranged from 3-12 ml depending on the size of the animal. Samples were kept chilled until they were able to be centrifuged at 3300 rpm for 10 min, 1-3 h post-collection. Separated plasma and spun red cells were aliquoted to cryovials and stored at −80 °C in a charged, dry nitrogen shipper (MVE Series DOBLE-20, New Prague, MN) until analyzed. Plasma (0.6 ml) and red blood cell (0.6 ml) samples were immediately sent to the Center for Applied Isotope Studies at the University of Georgia (Athens, GA) for bulk carbon and nitrogen stable isotope analysis, and plasma (1 ml) was sent to Nestle Purina Analytical Labs (St. Louis, MO) for fatty acid profile testing.
Stable Isotope Analysis. Blood samples (plasma and red blood cells) were analyzed for bulk nitrogen and carbon elemental compositions and stable isotope ratios of samples were determined using automated continuous-flow isotope ratio mass spectrometry 94 by the Center for Applied Isotope Studies at the University of Georgia. Lipid extraction and urea removal was performed on blood samples due to the potentially conflicting role of lipids and urea in isotope analyses of elasmobranch tissues 71 . For lipid extraction, subsamples (1000 µl plasma, 500 µl red blood cells) were treated with petroleum ether (5 ml), sonicated for 15 min, and the used solution was decanted 95 . This cycle was repeated for two additional rinses until the extracted material volume remained stable. Samples were then exposed to three cycles of rinsing with ultrapure deionized water (5 ml), sonication for 15 min, and decanting for the removal of urea 88 . Sample extracts were lyophilized overnight, homogenized, weighed, and then were combusted in tin capsules using an elemental analyzer (Model NA-1500, Carlo Erba, Milan, Italy) and injected into an isotope ratio mass spectrometer (Model Delta V, Thermo Scientific, Waltham, MA). Delta (δ) values were calculated based on the following equations where Vienna Pee Dee Belemite carbonate [V-PDB] and atmospheric nitrogen [Air-N2] are standards for carbon and nitrogen, respectively: and To normalize measured isotopes to an international scale, working standards (spinach and bovine liver) were regularly analyzed across a range of C and N masses and periodically checked in-house against primary international standards (e.g., National Institute of Standards and Technology). Routine measurements were precise to within 0.15‰ for δ 13 C and δ 15 N. Data are presented as mean ± SD. Fatty acid analysis. Plasma samples were sent to Nestle Purina Analytical Labs (St. Louis, MO) for FA profile analysis. Briefly, lipids were extracted using the Folch method 96 , which uses a 2:1 chloroform-methanol mixture. Fatty acids were esterified with methanolic sulfuric acid, taken up in heptane, and injected on a gas chromatograph (Model 6890, Aligent Technologies, Palo Alto, CA) equipped with a non-polar Equity-1 fused silica capillary column (100 m × 0.25 mm i.d., 0.2 µm film thickness), a flame ionization detector, and a split/ splitless injector. Helium was the carrier gas. The concentrations of individual FAs were calculated using an external standard method based on analysis of standards with known concentrations. The percentage for each FA was converted from the area of chromatogram peaks. Gas chromatograph results are typically subject to an error of up to ± 5% of individual component area. All FAs are expressed as percentage of total FA. Fatty acids are denoted as X:YnZ nomenclature, where X is the number of carbon atoms, Y is the number of double bonds in the carbon chain and nZ is the position of the first double bond from the terminal methyl end of the molecule. Cis or trans configuration of the double bonds are notated with a C or T, respectively. Data analysis. Data were evaluated for skewness, kurtosis, and normality using Kolmogorov-Smirnov test.
Data for δ 15 N and δ 13 C were not normally distributed and non-parametric Kruskal-Wallis tests were chosen for bivariate analyses along with Mann-Whitney tests to identify significantly different pair-wise comparisons. δ 15 N and δ 13 C data were power-and log-transformed, respectively, to normalize the distributions for each. Maximum likelihood linear mixed models were used to explore differences in nitrogen and carbon with sex, sampling time point, and an interaction of the two in order to determine whether sex differences varied across the sampling periods. Pairwise comparisons between the combinations of sampling time and sex were adjusted for multiple comparisons using the Šidák correction. A separate model was run for each tissue type and included a random www.nature.com/scientificreports/ intercept for individual animal to account for the between-subject heterogeneity. All testing was performed using Stata/SE statistical software (StataCorp, Stata Statistical Software: Release 14, College Station, TX). Multivariate non-parametric procedures were performed to describe differences in FA profiles based on 37 FAs that were present in proportions > 0.1% of the total. FA percentages were transformed using a square root transformation to avoid over-emphasis of extreme values prior to statistical analysis. Resemblance matrices were calculated using Bray-Curtis similarity measured between samples. We measured differences in profiles using a two-factor permutational multivariate analysis of variance (PERMANOVA) with sampling time point and sex as fixed factors. Similarity percentage routines (SIMPER) were employed to identify the FAs that contributed most to dissimilarities between males and females at each sampling point. Canonical analysis of principal coordinates (CAP) was then used to further discriminate FAs between a priori groups (males and females within each sampling time point), with FAs identified by SIMPER overlaid on the CAP ordination plot. Multivariate analyses were performed using PRIMER/PERMANOVA + software (version 7, Primer-E Ltd, Aukland, New Zealand). License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.