Fine-scale movement and habitat use of whitespotted eagle rays Aetobatus narinari in the Indian River Lagoon, Florida, USA

Despite being crucial to the conservation of batoids (skates and rays), assessments of fine-scale movements and habitat use of these taxa are lacking in the scientific literature. Here we used active acoustic telemetry to characterize habitat use and movement behavior of the stateprotected whitespotted eagle ray Aetobatus narinari in the Indian River Lagoon, Florida, USA. Seven mature individuals (6 males and 1 female) were individually tracked for a total of 119.6 h. Brownian bridge movement models of ray distribution showed the importance of habitats with high anthropogenic activity (i.e. boat traffic) such as inlets and channels, as well as clam aquaculture lease sites close to shore. This was supported by the significantly lower rates of movement in these habitats relative to other regions (offshore, open lagoon). Rate of movement significantly increased with temperature, suggesting that rays are more active during warmer periods. No tidal patterns in ray habitat use or distribution were evident. On average, rays used the deeper portions of the lagoon during the day and shallower portions during the night. While more extensive tracking is required to elucidate long-term movement patterns, this study is the first to characterize fine-scale habitat use by the whitespotted eagle ray in Florida while also identifying areas of potential interactions between this species and multiple anthropogenic threats.


INTRODUCTION
Estuaries and lagoons support a wide diversity of fishes and invertebrates due to high primary and secondary productivity (Beck et al. 2001). These complex systems are dynamic, varying considerably in abiotic factors (e.g. temperature, salinity, tidal range, sedimentation) that influence the diversity and behavior of animals in the system (Kennish 2002). Acquiring information on how large predators use these dynamic habitats is challenging, as many of these species are highly mobile; together, these factors make extended observations difficult. Fortunately, biotelemetry (i.e. fitting individuals with electronic tags and tracking movements) has emerged as a valuable tool to acquire habitat use data of highly mobile organisms in coastal ecosystems and has had demonstrated utility in animal conservation (Blumenthal et al. 2006, Schofield et al. 2007, Hays et al. 2016, Brooks et al. 2019. This approach can also help assess how anthropogenic disturbances common to coastal regions (e.g. habitat alteration, pollution, fishing, and boat traffic) can alter the structure and dynamics of specific populations and potentially reduce the sustainability of ecosystems within a coastal area (Kennish 2002, Curtis et al. 2013. In addition to understanding how an organism uses an ecosystem, comprehending the drivers of movement can also be vital for effective species conservation (Ogburn et al. 2017). Fishes in freshwater and marine ecosystems often move in response to a variety of biotic and abiotic factors (Jackson et al. 2001, Ortega et al. 2009, Cartamil et al. 2010, Schlaff et al. 2014. For example, the fine-and large-scale movements of many shark and ray species are attributed to abiotic factors (temperature, salinity, pH, tidal cycle), with the animals seeking out and/or remaining within the preferred environmental conditions (Schlaff et al. 2014). In particular, temperature is known to affect the physiology of ectothermic fishes, such as most elasmobranchs, by impacting functions such as digestion, reproduction, and growth (Bernal et al. 2012). Fine-scale movements to locate a preferred temperature range or behaviorally thermoregulating may be important for foraging in elasmobranchs (Matern et al. 2000, Sims et al. 2006, Hight & Lowe 2007, Speed et al. 2012. Moreover, understanding how abiotic factors influence the movement of a particular species can help model activity and distribution, as has been done with delta smelt Hypomesus transpacificus in the San Francisco Estuary, USA (Feyrer et al. 2011), deep sea fishes in the Azores (Parra et al. 2017), fishes in the Gulf of Mexico (Miller et al. 2018), and bull sharks Carcharhinus leucas in the Florida Everglades, USA (Heithaus et al. 2009). Modeling has also been used to understand how a species will respond to increased anthropogenic activities (Ivanova et al. 2020), which can help to create or improve conservation plans. However, these approaches require empirical data on habitat preferences in order to project a species' distribution and response.
Habitat use, fine-scale movement, and drivers of movement are also influenced by anthropogenic activities within an ecosystem. Understanding the spatio-temporal habitat use within highly modified coastal habitats (e.g. inlets, channels, docks, etc.) is essential for effective conservation, particularly in estuaries where anthropogenic activities have the potential to negatively affect populations (Kennish 2002, Waltham & Connolly 2011. For example, Werry et al. (2012) found that neonate and juvenile bull sharks on Australia's Gold Coast preferred natural to artificial habitats and argued that the decline in bull shark populations may be due to the vulnerability of this species to habitat modification associated with urbanization. Jennings et al. (2008) hypothesized that decreased juvenile lemon shark Negaprion brevirostris survival and reduced seagrass cover in Bimini, the Bahamas, was in part caused by dredging activity in the area. Anthropogenic alterations to ecosystems may also create more favorable habitats, altering natural movement and distribution patterns. For example, Barash et al. (2018) found several shark species aggregating more frequently near coastal power plant outflows in the Mediterranean Sea, likely due to warmer waters, when compared to other habitats. These alterations or preferences for anthropogenic habitats result in changes to behavior, habitat use, and movement that could lead to increased exposure to fishing and pollutants (Kennish 2002). Although some studies have documented the potential ramifications of using altered and anthropogenic habitats by top predators, data remain lacking for mesopredatory elasmobranchs such as rays.
The whitespotted eagle ray Aetobatus narinari, a mesopredatory batoid in the subtropical western Atlantic Ocean, has been documented in estuaries and lagoons throughout Florida (Gilmore 1977, Snelson & Williams 1981, Gilmore 1995, Bassos-Hull et al. 2014. Whitespotted eagle rays are listed as 'Near Threatened' on the IUCN Red List of Threatened Species, with a decreasing population trend due to the decline of mature individuals and fragmented populations (Kyne et al. 2006). While whitespotted eagle rays are protected from harvest throughout Florida, there have been no assessments of how the species interacts with anthropogenically derived features associated with developed coastal areas in the state. In Florida's Indian River Lagoon (IRL), examples of such features include channels (dredged inlet passes, intracoastal waterways), modified shorelines (docks, jetties, etc.), and shellfish culture beds. These features alter the natural benthos and create artificial habitats; however, the use of these artificial habitats by whitespotted eagle rays is not known.
The IRL has historically supported considerable shellfish fisheries, including a hard clam Mercenaria mercenaria fishery (FFWCC 2017). With declines in productivity over recent decades, much of the hard clam fishery is now dependent on aquaculture farms to produce seed, which is then planted on the natural bottom of the IRL for growout to market size. These growout sites are expansive in the Sebastian area of the IRL, yet the means by which potential predators such as whitespotted eagle rays interact with these modified habitats is unknown (Arnold et al. 2002). Collectively, clam lease holders have reported considerable interactions with whitespotted eagle rays and other predatory species over the years (E. Mangano pers. comm.), which has been further supported by presence data from passive acoustic receivers moored in select lease sites (B. C. DeGroot unpubl. data); however, the coarse scale of these data does not currently provide sufficient resolution to explain the nature of these interactions (i.e. traveling versus foraging). Understanding the interactions between whitespotted eagle rays and clamming activities in Sebastian will help sustain the last re maining components of this hard clam fishery and identify areas for successful bivalve restoration, which is a major goal of the IRL Council (IRLNEP 2018). Should these anthropogenic alterations create more favorable habitats for the species (e.g. enhanced foraging opportunities), interactions with humans may be increased and subsequently elevate the vulnerability of whitespotted eagle rays. Thus, an understanding of how this species uses natural and altered habitats is key to its population success in these coastal waters. Given that these habitats vary on fine spatial scales, high-resolution tracking approaches are required to obtain these important data.
The primary goal of this study was to employ active acoustic telemetry techniques to characterize whitespotted eagle ray habitat use, movement patterns, and interactions with anthropogenic habitats in the IRL. Although the fine-scale movement ecology of whitespotted eagle rays has been studied in Bermuda and the Bahamas using active acoustic telemetry (Silliman & Gruber 1999, Ajemian et al. 2012, neither study quantified movement trajectories or rates, nor did they examine habitat use in relation to anthropogenic activities. Additional goals of this study included using the associated habitat data collected during tracking to identify potential abiotic drivers of movement. Since previous studies found that whitespotted eagle rays exhibited diel movement patterns and evidence of tidal influences (Silliman & Gruber 1999, Ajemian et al. 2012, we hypothesized that habitat use and rate of movement would not be homogeneous throughout the day or tidal cycle.

Animal capture and tagging
Whitespotted eagle rays were collected in the IRL using 3 different nets: a 200 m long × 3 m deep (one 100 m panel of 15.24 cm stretch monofilament and one 100 m panel of 20.32 cm stretch monofilament) gillnet (N = 1), a nylon seine net (500 m long × 4 m deep) (N = 1), and a 35.56 cm stretch braided nylon tangle net (200 m long × 3 m deep) (N = 5). In all sets, rays were visually sighted from a vessel, after which the net was rapidly set around the animal in a circle of approximately 50 m diameter. The net was then en closed on the animal to facilitate capture. Once the animal was brought onboard, a series of length measurements (disc width [DW] and total length [from snout to trailing edge of anal fin]; cm) and the animal's weight (kg) were collected following Bassos-Hull et al. (2014). Individuals were tagged with an external rototag (Floy Tag), which was attached to the animal's dorsal fin. Additionally, 2 acoustic telemetry tags were surgically implanted into each ray's peritoneal cavity. Surgical instruments were kept in sterile packaging until washed with a Benz-All concentrate solution prior to the surgical procedure. Rays were then placed on their dorsal side, which induced tonic immobility as described by Henningsen (1994). The incision site for the acoustic transmitter was injected with a 2% lidocaine hydrochloride solution and swabbed with a 10% povidoneiodine swabstick. To avoid potential damage to reproductive organs located on the animal's left side, a 2−3 cm incision was made to the right of the animal's midline in the peritoneal cavity. A coded acoustic transmitter (V13-1H, 13 × 36 mm, 11 g in air, nominal transmission rate of 60 s) used for an additional study and a continuous acoustic transmitter (V13-1H, 13 × 36 mm, 11 g in air, 1 s transmission rate; Vemco) were inserted so the tags lay lateral to the internal organs. V13 continuous transmitters pinged on a unique frequency (60,63,75,78,81,84 kHz), and had a battery life of 23 d. Incisions were closed with braided absorbable sutures coated with polyglactin 910 (OS-6, 36 mm 1/2C, Reverse Cutting coated VICRYL Plus antibacterial needle) using 2 interrupted sutures, each secured with a surgeon's knot. After surgery, rays were monitored in an onboard holding tank to ensure proper recovery (spiracle respiration) before release (usually a maximum of 5 min). Handling and tagging procedures were compliant with local collection permits (FWC SAL-16-1785-SRP, FPS 07241710A, US FWS 41572-2017-07) and university-based IACUC protocols (FAU Animal Use Protocol A16-16).

Active tracking
Active tracking occurred in the IRL north of Sebastian Inlet (27.89°N, 80.49°W) to Wabasso, south of Sebastian Inlet (27.54°N, 80.29°W), including the nearshore waters of the Atlantic Ocean (Fig. 1). Seven individuals (denoted as R1−R7) were captured, tagged, and released in the IRL. Individuals were tracked continuously for up to 29 h with a VH110 directional hydrophone connected to a VR100 receiver and deck box (Vemco). Additionally, some animals were located between 3 and 20 d after the initial track and tracked again for up to 8 h; these secondary tracks are hereafter referred to as a 'relocation.' A custom-built aluminum pole mount secured the hydrophone to the tracking vessel (25' [ca. 7.6 m] flatbottom skiff) and facilitated tracking while underway (<5 knots). The VR100 receiver recorded latitude and longitude, decibel strength (dB), gain, and acoustic tag identification number every second when the animal was in range. Additionally, bearing (°) was recorded at 5 min intervals from an onboard GPS (Garmin ® GPSmap 78sc hand-held GPS). The position of the boat was assumed to be the position of the ray. The range and signal strength of the transmitters varied with depth, turbidity, and ambient noise. For each ray's initial track, the average gain was +4.0 and average dB was 73.3. Based on range testing (see Text S1 at www. int-res. com/ articles/ suppl/ n042 p109 _ supp.pdf), we assumed the tracking vessel averaged 19 m from the ray (Fig. S1). During tracking, environmental parameters including temperature (°C), dissolved oxygen (mg l −1 and %), salinity (ppt), and pH were collected from surface waters every 1−3 h for R1−R2 and every hour for R3−R7 using a YSI ® Professional Plus (Xylem).

Analysis
Raw VR100 detection data were filtered to encompass 10 ms before and after the programmed period of the tag. All spatial visualization and processing were conducted using ArcGIS desktop 10.3 (ESRI) with the Movement Ecology Tools for ArcGIS (ArcMET; Wall 2014). Within ArcMET, data were downsampled by averaging positions over 5 min intervals to better represent ray movement rather than the fine-scale boat maneuvering required to closely track the animal (Speed et al. 2013).
Rate of movement (ROM) was calculated from distances between consecutive 5 min positions. The average ROM by 5 min intervals was plotted by time of day to examine variability in ROM throughout the 24 h period, and a Scheirer-Ray-Hare test (extension of a Kruskal-Wallis test) was conducted to determine if ROM differed by time of day within habitat type. Additionally, to examine the potential for capture or tagging to affect ROM, a 1-way repeated measures ANOVA was conducted on the average 30 min ROM for the initial 5 h for all individuals in the study apart from R1. To examine potential effects of environmental conditions on whitespotted eagle ray movement, a generalized additive mixed model (GAMM) was built to determine if ROM differed with predictor variables: hour of the day, depth of the water column, salinity, temperature, tidal state, and habitat. GAMM analysis was conducted using the 'mgcv' library in R version 3.6.1 (R Core Team 2014). The presence of interactions between predictor variables was assessed during data exploration using co-plots. Temperature and salinity were measured hourly and assumed to remain constant for that hour. Tidal state for Sebastian Inlet (station ID 872204) was grouped into flood, high, ebb, and low tide with an hour for high and low tide based on NOAA's tides and current webpage (https:// tidesandcurrents. noaa. gov/ noaatidepredictions. html).
Habitat type was split based on both depth and benthic characteristics. Habitat type included the lagoon, the ocean, and channels such as the Sebastian Inlet. The channels had relatively deep water column depth inside the lagoon (>3 m) and the benthos was likely mud or sand. All areas within the IRL apart from channels were considered 'lagoon' habitat, which was shallow (<2 m), with benthic composition ranging from muck to sand and seagrass. Ocean habitat consisted of sandy benthos and hard-bottom habitats with depths ranging from 1−10 m. Both initial and subsequent tracks for each individual were included in the analysis for all rays apart from R1 due to short track duration (1 h). Collinearity between covariates was assessed using generalized varianceinflation factor (GVIF) scores; covariates yielding values higher than 3 were removed and GVIF scores were recalculated (Zuur & Ieno 2012). Track ID was added to the model as a random effect to avoid pseudo-replication and account for individual variation (Fig. S2). An auto-correlation plot was used to assess if there was serial correlation between residuals, where a value at time t is a linear function of the value at t − 1 (Zuur et al. 2009). Several distributions were considered for the response variable, including Gaussian, gamma, and log-normal. A circular smoother was used for hour of day. Smoothing splines were automatically optimized using crossvalidation in the 'mgcv' package. Second-order Akaike's information criterion corrected for small sample size (AICc) scores were used for model selection. The model with the lowest AICc score and most parsimonious model was selected (Burnham & Anderson 2003). Diagnostic plots including Q-Q plots, histograms of residuals, response versus fitted values, and linear predictors versus residuals were used to validate the chosen model.
For each track, a bottom depth was recorded every 5 min using a GPSMAP® 7612xsv sonar (Garmin) mounted to the hull of the tracking vessel. The bottom depths underlying all whitespotted eagle ray positions were converted to a raster file within ArcGIS. The average depth in each grid cell (0.179 km 2 ) was extracted, and the proportion of available habitats (all habitats used by all tracked rays) was calculated based on the following bins (<1.0, 1.0−1.9, 2.0−2.9, 3.0−3.9, 4.0−4.9, and >5.0 m). A selection index for underlying (or water column) depth was created based on the total available habitats used by all tagged rays throughout the tracking period. It was assumed that all habitats used were available to all tracked rays. A selection index (W i ) for each depth was calculated using the formula: (1) where O i = proportion of locations at depth i, and π i = proportion of available depths from all tracked rays. A selection index value >1 indicates selection for that depth and a value of <1 indicates avoidance (Manly et al. 2002). The habitat selection index was calculated in the R package 'adehabitatHS' with the function 'widesI,' which calculates a selection index and related chi-squared tests to examine habitat selection (Calenge 2006).
To better understand habitat use in relation to the channels/inlet, the lagoon, and the ocean, separate polygons outlining these features were created manually in ArcGIS. The proportion of downsampled positions overlapping each polygon was calculated for each individual for each portion of the day (day, night, crepuscular). Since tracking took place throughout the spring and summer (April−August), sunset and sunrise differed as much as 34 min. For these diel analyses, 'day' was defined as sunrise through evening civil twilight, 'night' as the period between evening civil twilight and morning civil twilight, and 'crepuscular' as morning civil twilight to sunrise and sunset to evening civil twilight. The proportion of positions in the channels/inlet, lagoon proper, and ocean out of the total positions was examined based on time of day.
To quantify the utilization distribution (UD) for each animal's track, the ArcMET toolset was used to calculate Brownian bridge movement models (BBMMs) for every initial track based on 5 min downsampled positions. The BBMM is a temporally explicit method for estimating the UD based on the probability of an individual being at a location along the movement path conditional to the distance and elapsed time between successive locations (Horne et al. 2007, Sawyer et al. 2009, Jachowski et al. 2013. The BBMM was chosen over more traditional kernel density estimates (KDEs) because BBMMs account for temporal autocorrelation and large spatial data sets, assumptions which KDEs violate (Horne et al. 2007). In addition, BBMMs account for spatial uncertainty through the Brownian motion variance or error associated with the estimated position of a tracked animal and irregular sampled data (Horne et al. 2007, Papastamatiou et al. 2013. The BBMM parameters were set as follows: maximum tolerated gap of 5 min, the algorithm automatically estimated the animal mobility variance using the 'Golden Section Search in One Dimension' routine, telemetry error based on the average gain and dB for each track, and 1 s integration time-step (Horne et al. 2007, Wall 2014. Core activity was defined by ≤50% UD contour while the 95% UD contour was used to define the total activity space.

Tagging
Seven adult whitespotted eagle rays (mean ± SD = 136.2 ± 10.6 cm DW) were manually tracked in the IRL between June 2017 and August 2018 (Fig. 1). All rays were tracked between the months of April and August (Table 1). The total duration of tracks per individual ranged from 1.0 to 31.2 h (17.1 ± 11.8 h, N = 7), yielding a total of 119.6 h of tracking and 1497 downsampled positions for all rays combined (Table 1). On 2 occasions, direct observations were made of tagged rays actively foraging, as described by Silliman & Gruber (1999). The first observation was of R3 on 20 April 2018 during slack low and rising tide from 16:34−18:24 h. The second observation was of R6 on 16 August 2018 during rising tide from 18:54−19:13 h. In both instances, the rays were observed foraging in a small area for prey at an average speed of 0.4 km h −1 before increasing speed to an average ROM of 1.3 km h −1 to transition to another location and repeat the process.

ROM and depth
The first 30 min of each ray's track was omitted from analysis based on the results of a 1-way ANOVA, which found significant differences in ROM by halfhour increment (F = 4.0621, p < 0.001). ROM was significantly higher in the first 30 min compared to all other intervals during the first 5 h of tracking (Fig. 2). The mean ± SD ROM among all whitespotted eagle rays was 1.4 ± 0.9 km h −1 . The average daytime, nighttime, and crepuscular ROMs for all tracks were 1.3 ± 0.9, 1.6 ± 0.9, and 1.5 ± 1.0 km h −1 , respectively (Fig. 3). The average 5 min ROM was 1.14 ± 0.84 km h −1 in channels, 1.56 ± 0.78 km h −1 in the ocean, and 1.52 ± 1.01 km h −1 in the lagoon. Average ROM varied significantly among habitats as determined by a Scheirer-Ray-Hare test (H = 46.47, p < 0.001). The average ROM also varied significantly between time of day (H = 28.11, p < 0.001), and there was a significant interaction with habitat and time of day (H = 15.88, p < 0.001). A Dunn's test also revealed significant differences when comparing the ROM between the lagoon and channels (p < 0.001) and between channels and the ocean (p < 0.001) but not between the lagoon and the ocean (p = 0.98). A post hoc Dunn's test revealed a significant difference when comparing ROM between day and night (p < 0.001) but not between day and crepuscular (p = 0.37) or crepuscular and night (p = 0.84). The interaction between habitat and time of day was primarily driven by the diel pattern in the channels. ROM in the channel habitat was slower during the day than at night (Fig. 4).
Four of the 7 whitespotted eagle rays spent over 50% of the tracking period in channel habitats. The most available habitat, by proportion, occurred in 1−1.9 m bottom depth. However, the most disproportionately used habitat was water with depths of 3−3.9 m (Fig. 5A). The proportion of time spent in these different regions varied with time of day. At night, a larger proportion of time was spent at depths <2 m, whereas during the day, more time was spent at depths >2 m (Fig. 5B) the overall selection index was significant (χ 2 = 133.43, p < 0.001). A selection index for depth by all tracked rays indicated significant preferences for the very shallow (<1 m, p < 0.001) and deeper (2−4.9 m, p < 0.001) portions of the lagoon (Fig. 6). Depths between 1 and 2 m were significantly avoided (p < 0.001) even though these depths were the most available habitats. The selection index for depths >5 m was not significant (p = 0.185). Environmental variables fluctuated within and among individual tracks (Table S1). GVIF scores in -dicated collinearity between temperature and salinity (GVIF >3). Temperature was selected over salinity and GVIF scores were recalculated, yielding new values that were <3 for all covariates. No interactions were evident in data exploration co-plots and were therefore not in cluded in the GAMM analysis. A log-normal distribution was selected based on the model diagnostic plots. The auto-correlation plot for initial model fits exposed a decline in serial correlation with increased lag. As such, an auto-regressive model of order 1 (AR-1) was incorporated in the model.
The model with the lowest AICc score included habitat type and temperature (Table 2) which were both statistically significant within the model (Table 3). ROM increased with temperature and the highest ROMs occurred offshore. By contrast, the animals moved most slowly in the channels (Fig. 7). ROM was not significantly affected by tidal state.

Utilization distributions
Overall, the inlet/channels had the highest proportion of detections during the day, and lagoon use increased during the crepuscular and nighttime portion of the day (Fig. 8A). Although no clear patterns were observed in whitespotted eagle ray habitat use with regard to tidal state, use of channel habitat was the highest during ebb and low tides (Fig. 8B). The BBMMs showed individual variability in distribution and habitat use. The BBMM for R1 was not included in analysis due to the short (1 h) track duration. The BBMM for R2 and R3 highlighted core areas of activity that occurred in and directly surrounding the channels of the intracoastal waterway (Fig. 9A,B). The core UDs for individuals R4 and R7 occurred in Sebastian Inlet (Fig. 9C,F), while those for R5 and R6 were centered in shallow waters parallel to land and adjacent to the channels (Fig. 9D,E). Although individual variation occurred, overlaying total activity space (95% UD) of all rays revealed overlap centered around the inlet and channel (Fig. 9G). Overlaying all individual BBMMs highlighted high re-use (i.e. top 50%) within inlets and channels or lagoonal waters adjacent to and parallel to those features (Fig. 9H). The sole female tagged in the study, R5, was the only individual to show habitat overlap with the hard clam aquaculture leases in the area. Percentage re-use by this female (R5) was relatively low al though re-use was as high as 50%.

DISCUSSION
While general data on the occurrence and habitat use of whitespotted eagle rays in the IRL were provided by Gilmore (1977Gilmore ( , 1995 and Snelson & Williams (1981), the present study is the first to delineate fine-scale habitat use in Florida waters and assess how abiotic factors influence ROM for the species. Additionally, with extended multi-day tracks, this study also provided a comprehensive and continuous fine-scale dataset for the species. These prolonged observations revealed affinities for habitat of considerable recreational and commercial importance.
Tracking individuals immediately upon release and then relocating individuals up to several days afterwards provided an opportunity to assess changes in whitespotted eagle ray movement behavior over time. Recovery times have been quantified for several elasmobranchs, but show significant variation among species and capture techniques. For example, juvenile lemon sharks caught on rod and reel had a maximum recovery rate of 0.58 h (Bullock et al. 2015), while Atlantic sharpnose sharks Rhizoprionodon terraenovae exhibited normal swimming patterns 1.5 h post release (Gurshin & Szedlmayer 2004). Accelerometers showed that the average recovery of adult blacktip sharks Carcharhinus limbatus caught on rod and reel was 10.54 ± 3.78 h (Whitney et al. 2016). The recovery rates for batoid fishes are less well known. A continu- ous acoustic tracking study of Hawaiian stingrays Hypanus lata showed increased ROM for the first 2 h of the track, assumed to be due to handling and tagging stress (Cartamil et al. 2003). Similarly, the results herein showed faster, sustained ROM immediately post re lease (as high as 9.0 km h −1 ), potentially indicative of capture, handling, and/or tagging stress. However, ROM normalized 0.5 h after release, suggesting that white spotted eagle rays recovered to normal swimming speeds during this time period. This information is useful to future studies that may incorporate devices such as acceleration data loggers to white spotted eagle rays to capture a suite of animal behaviors (Meyer 2017), and for post-release as sessment studies that analyze behavioral impacts from capture fisheries (Whitney et al. 2016). The active tracking approach used here permitted quantification of more probable movement rates than those previously attained in other studies. The mean ± SD ROM among all whitespotted eagle rays in this study was 1.4 ± 0.4 km h −1 , substantially higher than the rates identified by Ajemian & Powers (2014), who documented maximum ROM of 0.56 ± 0.08 km h −1 for males and 0.28 ± 0.1 km h −1 for females, measured from towed-float satellite telemetry. Consideration of tag type may account for this difference. Position estimates from satellite telemetry require the animal to be close to the surface, and position solutions (via Doppler shift algorithms) are often associated with higher errors (at best 400 m) (Douglas et al. 2012). Thus, satellite telemetry may not capture the full extent of whitespotted eagle ray movements, especially finerscale transitions. As such, it is likely that the ROM documented here more closely approximates the actual capabilities of this species, which can help calibrate dispersal models and finetune satellite geolocation tools that require a priori knowledge of speciesspecific movement rates.
As with other studies examining the role of abiotic factors in the movement of sharks and rays (Matern et al. 2000, Sims et al. 2006, Cartamil et al. 2010, we found that ROM was significantly affected by temperature. As previously mentioned, temperature is a well-documented abiotic factor that affects the physiology of ectotherms (Bernal et al. 2012, Schlaff et al. 2014, and it is  Table 2. Top 5 most suitable generalized additive mixed models for investigating the rate of movement (ROM) of whitespotted eagle rays in response to the chosen covariates. 'f': covariate treated as a factor; 's': covariate treated with a smoothing spline. All models included a random effect for track ID and an autocorrelation structure to account for serial correlation in the data. Bolded model indicates the optimal model. AICc: Akaike's information criterion corrected for small sample size  thought that many elasmobranch species actively seek out preferred temperatures throughout the day (Schlaff et al. 2014). For example, lesser spotted dogfish Scyliorhinus canicula have shown diel differences in depth use in Ireland, with more time spent at depth in cooler waters during the day, whereas shallow warmer waters were occupied at night (Sims et al. 2006). In the current study, ROM steadily increased with rising temperature, which may be indicative of movement to locate a preferred temperature range (Matern et al. 2000, Cartamil et al. 2010), behavioral thermoregulation (McLaren 1963, Matern et al. 2000, predator avoidance (Heupel & Heuter 2002, Collins et al. 2007, Heithaus et al. 2009), or a combination of several of these factors. Furthermore, at warmer temperatures, more enzymes have reached their activation energies, allowing for more physical activity which may contribute to the increased ROM with increasing temperature (Halsey et al. 2015). Additional tracking studies that incorporate continuous transmitters with temperature and pressure sensors may help elucidate the role of the behaviors observed herein. Additionally, we did not find that ROM was significantly affected by tidal state as has been reported for other elasmobranch species (Schlaff et al. 2014). Furthermore, we observed no pattern with habitat use by tidal state. This may be attributed to the relatively low sample size, and future studies should continue to evaluate the effect of tidal state on whitespotted eagle ray movement. Whitespotted eagle rays are benthic predators, consuming infaunal and epifaunal bivalves and gastropods (Schluessel et al. 2010, Ajemian et al. 2012, Serrano-Flores et al. 2019. These benthic prey have been documented in the IRL (Hesselman et al. 1989, Arnold et al. 2000 and generally prefer shallow (<2 m) regions of lagoons that facilitate faster growth (Arnold et al. 1991). Direct observations of feeding behavior during this study (e.g. R3 and R6 foraging at dusk in the shallow flats <1 m) showed that whitespotted eagle rays increased speed to move to a new location, searched the area at a slower rate, and then again increased speed to transition to another location and repeat the process. Based on this information, foraging would be associated with an increased ROM relative to a resting ROM, but a decreased ROM relative to predator avoidance or transitioning between habitats. Tracking data revealed that whitespotted eagle rays move most rapidly in ocean habi- tats, where risk of predation may be elevated, while they moved at the slowest speeds in the channels/ inlets, and at intermediate speeds in the lagoon. Additionally, the ROM was higher at night, when rays were spending a larger proportion of time in the lagoonal habitats, than during the day. This increased use of shallower lagoonal habitats at night, where bivalves occur at higher densities, may suggest that whitespotted eagle rays are transitioning between and foraging in these habitats. Consuming prey in this shallow, warm water and resting in cool water, such as inlets (where the ROM for whitespotted eagle rays is lowest), has been hypothesized to minimize metabolic costs for related bat rays Myliobatis californica; Matern et al. 2000). Further supporting the findings of Matern et al. (2000), white spotted eagle rays disproportionately selected very shallow (<1 m) and deeper (2−4.9 m) areas of the lagoon. Variation in ROM among habitat types and proportion of time spent in different habitats may also be in re sponse to predator avoidance, as several shark species are more active at night (Lowe et al. 1996, Whitney et al. 2007) and likely access lagoonal waters via inlets and channels of the IRL (Gilmore 1977, Snelson & Williams 1981. However, understanding predator avoidance behaviors would require additional fine-scale habitat use assessments for large coastal sharks in the IRL. Additionally, whitespotted eagle ray behavior was not always discernable due to low water clarity and/or the ray's position in the water column. To help elicit behavioral patterns, future fine-scale movement studies should also consider implementing accelerometer data loggers to provide additional insight into behavior at these locations and depths. As determined by the BBMM, whitespotted eagle rays did not re-use the interior portion of the lagoon and may have used this area to transition between habitats. In contrast, both the BBMM and habitat selection indices highlighted the importance of highly modified coastal habitats to all whitespotted eagle rays. The BBMM revealed the importance of channels and inlets for whitespotted eagle rays, supporting Ajemian et al. (2012), who observed large groups aggregating around Flatts Inlet in Bermuda. Similarly, Silliman & Gruber (1999) also deemed deep (3−6 m) channels with high currents as core areas for whitespotted eagle ray resting and refuging. Understanding channel use is crucial to evaluating risks and potentially developing strategies to mitigate negative impacts to the species (Cloyed et al. 2019), as both channel and inlet habitats have high levels of anthropogenic activity (e.g. boating, fishing) and are prone to coastal development impacts from dredging (Kennish 2002, Waltham & Connolly 2011. Additionally, these high traffic areas experience increased noise and chemical pollution (Pirotta et al. 2019). More over, in the present study whitespotted eagle rays spent a larger proportion of time in the channel and inlet habitats during the day, when boating and fishing activity tends to be highest. Therefore, anthropo genic impacts associated with high vessel traffic pose a potential threat to this species. Whitespotted eagle rays entangled in fishing line and crab traps as well as rostra pierced by fishing hooks have been reported in the region (M. J. Ajemian pers. obs.). As coastal populations and development increase, there is more potential for whitespotted eagle rays to be negatively affected by these activities. Additionally, intense coastal development such as dredging, construction, and pollution have been linked to habitat alteration (Lotze et al. 2006, Stump 2013, Wenger et al. 2018, which may change the abundance and distribution of whitespotted eagle rays, as has been documented with shark species in degraded habitats (Jennings et al. 2008, Werry et al. 2012, Curtis et al. 2013, Barash et al. 2018. As whitespotted eagle rays already display an affinity for these modified habitats, increased interactions with humans and added pollution and/or disturbances could result in changes to the species' movement patterns and health. Ultimately, such human-induced habitat alterations could reduce the overall productivity of estuarine areas and, with time, exacerbate pressures already facing populations of aetobatid rays (e.g. bycatch, entanglement; Kyne et al. 2006).
UDs exhibited by whitespotted eagle rays can also help inform local bivalve aquaculture farms and restoration efforts. Productivity of local hard clam leases in east-central Florida has steadily declined over the years (Arnold et al. 2002, E. Mangano pers. comm.). Be cause durophagous stingrays are attracted to high-density bivalve shellfish habitats where they can have measurable effects (Peterson et al. 2001, Fegley et al. 2009, Ajemian et al. 2012, the habitat use data collected herein could help potentially mitigate negative interactions with grow out activities where clam densities far exceed ambient levels. For example, due to whitespotted eagle ray affinity for shoreline and channel habitats, the development of new leases or bivalve restoration pro grams might consider concentrating efforts towards the interior of the lagoon, where less area re-use was documented. However, long-term monitoring of predator distribution may be required given the temporal limitation of this study.

Limitations and future considerations
The individuals tracked in this study mainly consisted of adult males. As such, our understanding of whitespotted eagle ray habitat use in this system is likely biased, given that elasmobranchs oftentimes exhibit different habitat use patterns by sex. For ex ample, female elasmobranchs generally show higher residency in lagoon habitats where parturition sites are located (Wearmouth & Sims 2008, Schlaff et al. 2014). Thus, our study may underestimate the level of interactions occurring between whitespotted eagle rays and anthropogenic activities in the IRL due to the more transitory behavior of males. Additional tagging and tracking studies that encompass all sizes and sexes throughout all seasons are needed for a more comprehensive study.
Another limitation of this study was our understanding of depth use by the tagged individuals. Although there was minimal error in positional estimates collected during active tracking (±19 m), depth was extracted from the vessel's on-board depth sounder to gain an estimate of water column habitat use. However, future studies should consider the incorporation of pressure sensors into acoustic tags to ascertain the actual position of a ray in the water column. Such data could provide further information on where the animal may be interacting with the benthos (e.g. clam growout sites), as well as better ascertain collision risk to vessel traffic given the preference of this species for channel habitats. In both cases, we recommend that future researchers use high-resolution pressure sensors given the generally shallow habitat (mean = 1.9 m) of the region.
Tracking durations (<29 h) in our study were limited by the memory of the VR100 receiver deck box; however, relocation attempts were conducted to verify behavior and provide additional tracking data. Tracking was affected and limited by weather (e.g. wind, rain, and lighting) as well as local vessel traffic. Longer tag periods or an additional VR100 receiver deck box are recommended to facilitate prolonged tracking durations, including additional nighttime tracking. More extensive datasets may provide a greater understanding of diel differences in ROM and habitat use and reveal additional drivers of movement, particularly in relation to foraging be haviors. Although environmental data were taken throughout the tracking periods, more frequent sampling may improve future GAMM analyses. Lastly, to help confirm behavior, future studies should consider incorporating accelerometry sensors and/or action cameras into tag packages while remaining mindful of tag package to animal weight ratios.

Conclusions
Despite the limited sample size and tracking durations, consistent patterns emerged among the 7 tagged individuals. Whitespotted eagle rays of the IRL system exhibited diel patterns in distribution and habitat use. Individuals spent a larger proportion of time in the channels and inlet during the lighter and warmer portions of the day and used shallower depths during the cooler and darker portions of the day, supporting findings of Ajemian et al. (2012) in Bermuda. In general, whitespotted eagle rays had faster ROMs while in the ocean and lagoonal habitats and slower ROMs in channels and inlets. This information suggests that whitespotted eagle rays may spend more time foraging at night in the shallow waters of the lagoon than during the daytime. Because most clammers work on lease sites during the day, interactions between whitespotted eagle rays and growout sites may therefore be underestimated. Additionally, this study was the first to document high re-use of anthropogenic habitats by this species. Overall, this research increases our understanding of fine-scale movement patterns and UDs of whitespotted eagle rays, which will help inform statewide conservation plans for the species and provide critical information to hard clam aquaculture farmers and restoration managers for the successful production of bivalves in the area.