Does space use behavior relate to exploration in a species that is rapidly expanding its geographic range?
McCune KB1*
Ross C2
Folsom M2
Bergeron L1
Logan CJ2
23 September 2020
Affiliations:
- University of California Santa Barbara
- Max Planck Institute for Evolutionary Anthropology
*Corresponding author: KB McCune (kelseybmccune@gmail.com)
Cite as: McCune K, Ross C, Folsom M, Bergeron L, Logan CJ. 2020. Does space use behavior relate to exploration in a species that is rapidly expanding its geographic range? (http://corinalogan.com/Preregistrations/gspaceuse.html) In principle acceptance by PCI Ecology of the version on 23 Sep 2020 https://github.com/corinalogan/grackles/blob/abc201d3bb05daa31fa85225747fb9a0682a3e33/Files/Preregistrations/gspaceuse.Rmd.
This preregistration has been pre-study peer reviewed and received an In Principle Recommendation by:
Blandine Doligez (2020) Explore and move: a key to success in a changing world? Peer Community in Ecology, 100058. 10.24072/pci.ecology.100058
- Reviewers: Joe Nocera, Marion Nicolaus, and Laure Cauchard
ABSTRACT
Great-tailed grackles (Quiscalus mexicanus) are rapidly expanding their geographic range (Wehtje 2003). Range expansion could be facilitated by consistent behavioural differences between individuals on the range edge and those in other parts of the range (Duckworth and Badyaev 2007; Lindström et al. 2013). Movement behaviors in particular are thought to be important for expanding and invasive species (Bubb, Thom, and Lucas 2006). There is evidence for a relationship between individual exploratory traits and dispersal [the permanent movement an individual makes from its birth site to the place where it reproduces; Greenwood and Harvey (1982)], but it is still unknown whether individual differences in exploration relate to adult daily movement patterns within the home range (“space use”). Evidence suggests that daily space use behavior can vary consistently among individuals (Hertel et al. 2020), but no studies have examined the relationship between space use behavior and measures of exploration, or space use behavior in different populations along an expanding range. Here we will study grackles from 3 populations that span the extent of the current range - Central America (their original range), Arizona (middle of the northern expanding edge), and northern California (near the northern edge of their range). We will test whether performance on an exploration task in captivity relates to subsequent space use behavior in the wild measured using home range size, the autocorrelation of step length (distance between two sequential observations) and turning angle for each individual over time (Pacheco-Cobos et al. 2019), and the repeatability of each individual’s occurrence in particular geographic locations. Results will inform whether individual differences in space use behavior are associated with consistent individual differences in exploration, and whether daily space use patterns differ in populations at different points on an expanding range. If space use behavior correlates with measures of exploration, then space use data could be used to inform conservation management strategies (e.g. identify which individuals are likely to remain in new or restored habitat after a translocation (May, Page, and Fleming 2016)) in species where it is not logistically feasible to a priori measure exploration in captivity. If daily space use patterns differ in populations at different locations on the expanding range, then future studies could investigate whether space use and exploration are coupled with other traits such as long-distance natal dispersal to form an “invasion syndrome” (Chapple, Simmonds, and Wong 2012).
INTRODUCTION
Problematic range expansions of invasive species are occurring across the globe. However, not all species that invade novel areas outside of their traditional range become established (Chapple, Simmonds, and Wong 2012; Hayes and Barry 2008) and it is still unclear what characteristics facilitate successful invasions (Fogarty, Cote, and Sih 2011). These characteristics are important for predicting potential invasions, especially since invasive species are implicated as a leading cause of biodiversity loss (Butchart et al. 2010; Clavero and Garcı́a-Berthou 2005). One hypothesis is that the composition of behavioral types might facilitate geographic range expansion (Carere and Gherardi 2013) because certain individuals may possess traits such as exploration and aggressiveness that make them more likely to succeed when venturing into new habitats and competing with heterospecifics (Cote et al. 2010). Evidence is already accumulating that links consistent individual differences in behavioral types and invasion success. For example, Duckworth and Badyaev (2007) found more aggressive behavior in western bluebird males on the range edge, suggesting that range expansion was facilitated by aggressive males dispersing further and displacing less aggressive mountain bluebirds. Additionally, invasive cane toads on the invasion front in Australia spent longer in the dispersal mode and moved greater distances than toads tracked several years later in the same location (Lindström et al. 2013).
Within-species variation in the ability (movement) and motivation (exploratory tendency) to encounter conspecifics, novel foods, and novel food sources could be a limiting factor in successful species range expansions (Spiegel and Crofoot 2016). In novel areas, the occurrence of conspecifics, food, predators and other environmental factors may not be as easily detectable or recognizable, and may be distributed differently across the landscape than in core areas of the range. Although individuals with exploratory phenotypes may be more successful at colonizing new habitat, exploratory individuals are also at higher risk of predation (Stuber et al. 2013), and could be less likely to find local food sources (Overveld and Matthysen 2009). Consequently, for establishment in new areas, individuals that exhibit a range of exploratory behavior are needed, and the interaction between space use and exploratory tendency is likely important for finding novel foods and food sources. Additionally, while dispersal [the permanent movement an individual makes from its birth site to the place where it reproduces; Greenwood and Harvey (1982)] is necessary to initially invade the novel habitat, subsequent daily space use could determine establishment success. For example, on the range edge conspecific density might be lower and individuals may need to use space differently (Bubb, Thom, and Lucas 2006). However, current research on invasion success and movement behavior focuses heavily on dispersal, and while dispersal and exploratory tendencies have been shown to be associated (Cote et al. 2010), we do not know how exploratory tendency influences space use patterns in the daily lives of invading individuals.
Space use behavior is expected to be influenced by internal states like exploratory tendency and hunger, as well as the non-random distribution of available habitat and resources (Nathan et al. 2008). Space use can also consistently differ among individuals (Hertel et al. 2020), which indicates that each individual has distinct preferences for how, when, and where to move within its home range. Traditional analyses of animal space use required spatial and temporal independence of data points for statistical analysis (Swihart and Slade 1985), yet movement data are unlikely to meet these criteria. Spatial and temporal autocorrelation occurs when the position of an individual at a given time is tightly linked to its position both before and after, including cases where individuals are repeatedly found in the same locations across time. This autocorrelation is an intrinsic component of space use behavior and eliminating it can reduce biological relevance and obscure relationships with behavioral types (Dray, Royer-Carenzi, and Calenge 2010). Therefore, the autocorrelated nature of movement paths could be important to illuminate the relationship between individual differences in exploratory tendency and daily space use.
Great-tailed grackles (Quiscalus mexicanus, hereafter “grackles”) are rapidly expanding their geographic range (Wehtje 2003). They are invasive because they meet the criteria for the establishment and spread stages of invasion (Blackburn et al. 2011), and they are considered a pest in some areas where they reduce fruit crop yields (Glahn, Palacios, and Garrison 1997). The nature and level of ecological and social factors grackles experience may vary in importance between populations. Generally, this species is strongly associated with human-modified landscapes and is able to take advantage of a variety of human foods (e.g. crops, at our outdoor cafes, and out of our garbage cans) in addition to foraging on insects and on the ground for natural food items (K. Johnson and Peer 2001). Furthermore, they exhibit variation in territorial social behavior in both the breeding and non-breeding seasons. During the non-breeding season, this species forages in smaller groups and communally roosts in larger groups (K. Johnson and Peer 2001). During the breeding season, one or more males defend a territory where multiple females place their nests within that territory to raise the young (K. Johnson et al. 2000). Roaming males are also present and can obtain extra-pair copulations with females on other males’ territories (K. Johnson et al. 2000).
In this investigation, we aim to understand whether exploratory tendency, measured in captivity, is associated with space use behavior measured in the wild in grackles from three populations that span the current range: Central America (their original range), Arizona (middle of the northern expanding edge), and northern California (near the northern edge of their range). Exploration, measured here following the protocol described in McCune et al. 2019, is interpreted as an individual’s response to novelty, such as novel environments or novel objects (Réale et al. 2007), to gather information that does not satisfy immediate needs (Mettke-Hofmann, Winkler, and Leisler 2002). Previous studies on birds have found a relationship between some movement behaviors (home range size and natal dispersal) in the wild and exploration measured in captivity using novel environment (Minderman et al. 2010; Dingemanse et al. 2003) and novel object (Mettke-Hofmann, Winkler, and Leisler 2002) tests. Therefore, these measures have the potential to be relevant to grackles.
After measuring exploratory tendency in captivity, we will release grackles and use telemetry to measure their space use behavior. We will quantify a typical measure of space use that controls for autocorrelation (i.e. home range size), but we will also investigate the behavioral relevance of autocorrelation in space use with two relatively new methods. The first method will describe individual differences in path-level movement behavior by analyzing the autocorrelation of step length (distance between two sequential observations) and turning angle for each individual over time (Pacheco-Cobos et al. 2019; Hertel et al. 2020), while the second method will describe individual differences in spatial preferences by analyzing the repeatability of each individual’s occurrence in particular geographic locations.
With these data we will test two separate hypotheses. First, we will test whether grackles’ performance on exploration tasks in captivity is related to the space use metrics of the same individuals that we quantify from their movement behavior in the wild. Second, we will test the hypothesis that the space use behavior of wild grackles systematically varies with population location along the extent of the range. These results will have implications for invasive species management (Carere and Gherardi 2013), but also could be used to inform the conservation of threatened and endangered species [e.g. identify which individuals are likely to remain in new or restored habitat after a translocation; May, Page, and Fleming (2016)]. This will be especially relevant for species where it is not logistically feasible to a priori measure exploration in captivity.
A. STATE OF THE DATA
This preregistration uses secondary data: data that are already being collected for other purposes (GPS points in hypothesis 3 and home range sizes in prediction 3 in the flexibility and foraging preregistration). Originally we attached radio tags to grackles released from the aviaries to ensure that we could find their nest sites and track measures of foraging behavior and reproductive success for these individuals for which we have an extensive amount of behavioral data from aviary tests. Now we plan to additionally use the radio tags to collect data for this space use preregistration. This preregistration was written in June 2019, while at the same time increasing the number of GPS points taken per tracking session per bird to provide enough data for the analyses here, and submitted in September 2019 to PCI Ecology for pre-study peer review. Reviews were received in December 2019 and we revised and resubmitted in March 2020. A second round of reviews was received in July 2020 and we revised and resubmitted in August 2020. A third round of reviews was received on September 1, 2020 and we revised and resubmitted on September 15, 2020. A fourth round of reviews was received on September 17, 2020 and we revised and resubmitted on September 30, 2020.
B. HYPOTHESES
H2: Space use behavior will vary among grackles from our three study populations located along different points in the geographic range of this species (core, middle of expansion, and range edge; Table 1). These populations are theoretically connected, however actually moving between two of our field sites within a few grackle lifespans is unlikely due to the large distances between field sites and two geographic barriers (the Sierra Nevada and Sierra Madre mountain ranges, and the high elevation areas of Mexico).
Prediction 2: Home range size will increase, autocorrelation of step lengths and turning angles will decrease (i.e. grackle movement behavior will be less predictable), and grackles will use a greater variety of spatial locations as the geographic distance from the original center of the range increases. Specifically, the grackles sampled from our site on the edge of the range (northern California), will have larger overall home range sizes, exhibit more variety in step lengths and turning angles, and use a greater variety of spatial locations than the sample of grackles in the core of the range (Central America). Grackles in the sample from the middle of the expanding range will be intermediate in space use (Arizona). Such population differences in space use behavior may relate to range expansion because some of the individuals on the leading edge of the range may use more space and move longer distances (Duckworth and Badyaev 2007). However, larger-scale sampling of grackle groups across the strata of the expansion front and core range would be needed to more robustly validate the hypothesis that our cross-site differences are indicative of a broader pattern driven by the location of the expansion front.
Prediction 2 alternative 1: Grackles sampled on the edge of the range will have smaller overall home range sizes, high autocorrelation in step length and turning angle (i.e. movement behavior will be more predictable), and consistently use the same spatial locations compared to grackles sampled in the middle or core of the current range. This would suggest that suitable habitat may be distributed in small patches, that novel habitats at the edge of the range may have high predation on grackles that use more space, and/or that individual grackles specialize on certain novel habitat types that are patchily distributed.
Prediction 2 alternative 2: We will find no difference across the geographic range in the space use behavior of the grackles sampled. This would suggest that, on average, all grackles may use the same amount of space, or that there is a similar distribution of individual differences in space use in each population. Alternatively, it could indicate that we did not detect differences because we measured adults rather than juveniles. Grackles sampled in different populations may converge on similar space use behavior during development, or juvenile grackles may disperse further on the edge of the range. However, we are not able to detect these differences with our data, which is primarily from adults.
Table 1. Population characteristics for each of the three field sites. The number of generations at a site is based on a generation length of 5.6 years for this species (International 2018) and on the first year in which this species was reported to breed at the location (Wehtje (2003) for Arizona, Steve Hampton’s pers. comm. reported in Pandolfino, Deuel, and Young (2009) for Woodland, California). The first confirmed nest sighting in Woodland, California was reported in the Yolo Audubon Society’s newsletter The Burrowing Owl (July 2004), which Steve Hampton shared with Logan. For Central America, there is no data on the first year in which they started breeding because this species originates in this region, therefore we used the age of the species: 800,000 years (N. K. Johnson and Cicero 2004).
d <- read.csv(url("https://raw.githubusercontent.com/corinalogan/grackles/master/Files/Preregistrations/gspaceuse_Table1.csv"), header=F, sep=",", stringsAsFactors=F)
dat<-data.frame(d)
library(reactable)
reactable(dat, highlight=TRUE, bordered=FALSE, compact=TRUE, wrap=TRUE, resizable=TRUE,
columns = list(
V1 = colDef(name="Site"),
V2 = colDef(name="Range position"),
V3 = colDef(name="Breeding since"),
V4 = colDef(name="Number of years breeding"),
V5 = colDef(name="Average number of generations"),
V6 = colDef(name="Kilometers to nearest field site"),
V7 = colDef(name="Nearest field site"),
V8 = colDef(name="Citation")
))C. METHODS
Planned Sample
Great-tailed grackles are caught in the wild, given colored leg bands in unique combinations for individual identification, and released at their point of capture. The color-marked grackles in this study have one of two different backgrounds: those that do not have radio tags and those that do. First, we opportunistically track color-marked grackles that do not have a radio tag (and thus have not spent time in the aviaries) to compare whether time spent in the aviaries is related to space use behavior. When a color-marked bird is encountered, researchers track it for 20-90 minutes, recording the spatial location every one minute. If the bird goes out of view, researchers attempt to find it again for 15-30 minutes before moving on. Because these data are opportunistic, we do not attempt to balance for sex, but we aim to follow at least 20 non-tagged individuals in each population.
Second, we applied radio tags to all aviary-tested birds (estimated 20 individuals per population). These subjects are primarily adults who we attempt to balance for sex, and who spent up to six months in an aviary while they participated in behavioral choice tests (see Logan et al. 2019 for details) and individual differences assays, including measures of exploration in captivity (see McCune et al. 2019 for details), as part of other research projects by this lab. For details about the captive environment, please refer to the preregistration associated with this part of the research: McCune et al. 2019. Cognitive and behavioral traits are often not fully developed in hatch year birds (i.e. Zucca, Milos, and Vallortigara 2007). For our aviary cognitive test battery, we avoided the potentially confounding variable of stage of cognitive development by only taking known adult grackles into the aviaries. We identified grackles as adults using eye color, where hatch year birds have brown eyes, but second year and older birds have yellow eyes. While it is possible that cognitive traits continue to develop in second year and older birds, it is impossible to distinguish grackle age after the bird’s first year (K. Johnson and Peer 2001).
Before the aviary grackles are released, they are fitted with VHF radio tags (Lotek PipLL (model Ag386), Advanced Telemetry Systems (model A2455) or Holohil Systems Ltd. (model BD-2)) so we can track space use behavior using radio telemetry. Radio tags were initially attached to the grackles by gluing them to their backs (K. Johnson and Peer 2001; Mong and Sandercock 2007), however these did not stay on for very long. Therefore, we now use a leg loop harness (methods as in Rappole and Tipton 1991) made from sutures and secured with crimp beads and cyanoacrylate glue (Vicryl undyed 36in sutures, item number D9389 at eSutures.com; 0.5mm diameter, absorbable so they fall off after one to four months).
After release, an experimenter will find and follow (each session is called a “track”) each tagged grackle at least four times per week. On each of these tracks, the experimenter follows the focal individual for approximately 1.5 hours, recording a GPS point every one minute, regardless of whether the bird moved (Cushman, Chase, and Griffin 2005). Additionally, we aim to balance tracking data equally during morning and afternoon time periods for all grackles. We hoped to track the same individuals across the breeding and nonbreeding seasons. However, the leg loop harnesses often degraded and fell off after 1 - 4 months. Therefore, we have few data points on the same individuals in both seasons and we will instead compare space use in different individuals (rather than within individuals) across seasons. Researchers maintain a distance of at least 30 m and observe the bird with binoculars so the grackle’s behavior is not influenced or artificially changed. If the grackle alarm calls while oriented towards the researcher (indicating the researcher’s presence affected the grackle’s behavior), all tracking on that individual is stopped for the day. To ensure we capture all locations the individual visits and not just those where they are most easily seen and followed, tagged grackles that move out of sight during tracking are searched for with telemetry until they are found again.
Exploratory tendency is an intrinsic factor that could be related to space use behaviors. However, there are also likely alternative variables that may relate to space use behavior in wild grackles that we must control for by including them as covariates in our models. First, we measure energetic condition (described in Berens et al. 2019) to account for differences in the physiological mobility that may limit an individual’s space use behaviors (Nathan et al. 2008). Secondly, we measure habitat characteristics such as human food sources and available breeding habitat (described in Logan et al. 2019) because these factors of the external environment will affect where grackles choose to move or spend time (Nathan et al. 2008).
Conspecific density has also been shown to affect home range size in other bird species (Flockhart et al. 2016; Garabedian et al. 2018). To control for the possibility that home range size may vary among our populations due to conspecific density rather than exploratory traits, we will use point count surveys to measure grackle population density. We will place 225 point count stations across the landscape encompassing each population (Tempe, AZ; Woodland, CA; Gamboa, Panama). For each study population, the first central point will be randomly placed within 500 m of the center of the study area. The remaining points are placed in a 500 m grid pattern extending out from this central point. In total, the sample area will cover an area that is 7 km by 7 km. Each point will be visited once during the non-breeding season (Sep-Mar). During the survey, researchers will record all grackles visually and aurally detected for six minutes.
Sample size rationale
We test as many birds as we can during the approximate five years of this study given that the birds are only brought into the aviaries during the non-breeding season (approximately September through March). It is time intensive to conduct the aviary test battery (2-6 months per bird at the Arizona field site), therefore we approximate that the minimum sample size for captive subjects will be 60 across the three sites (approximately 20 birds per site with the aim that half of the grackles tested at each site are females). We catch grackles with a variety of methods (mist nets, walk-in traps, and bow nets), some of which decrease the likelihood of a selection bias for exploratory and bold individuals because grackles cannot see the traps (i.e. mist nets). Once released, we will primarily track the space use behavior of these ~60 grackles that have radio tags. As described above, we will also opportunistically collect GPS point locations non-tagged grackles to determine whether grackles that were previously in the aviary have different space use behavior from non-aviary-held grackles after their release. We will attempt to match the sample size of aviary birds, and in our Arizona population we currently have over 20 points [the minimum number for reliably calculating home range size; Noonan et al. (2019)] for 33 individuals that have never had radio tags. We aim to acquire more than 20 points on at least 20 non-tagged grackles in the other two populations as well. Additionally, we attach radio tags to birds that are released early because of their lack of willingness to participate in aviary tests (currently 5 individuals) to determine whether space use behavior differs between participatory and non-participatory grackles.
Data collection stopping rule
We will stop collecting GPS location data on tagged and non-tagged birds when home ranges are fully revealed for data collected in both breeding and non-breeding seasons. To determine at what point home ranges have been fully revealed, we will calculate the asymptotic convergence of home range area as in Leo et al. (2016). We will test home range asymptotic convergence for breeding season and non-breeding season movements separately (breeding season: Apr - Aug, non-breeding season: Sep - Mar).
Open materials
Protocols:
- Exploration protocol for exploration of new environments and objects, boldness, persistence, and motor diversity.
- Radio tracking protocol for attaching radio tags and collecting GPS points using radio telemetry.
- Point count protocol for measuring grackle population density in the study area.
Open data
When the study is complete, the data will be published in the Knowledge Network for Biocomplexity’s data repository.
Randomization and counterbalancing
There is no randomization in this investigation. The order of the exploration tasks is counterbalanced across birds (see the separate preregistration for details). The time of day that we collect GPS point locations is counterbalanced within and across birds to account for potential variation in movement behavior arising from daily circadian rhythms.
Blinding of conditions during analysis
No blinding is involved in this investigation.
Summary of methods for measuring exploration McCune et al. 2019
We adapted commonly used methods to test exploratory tendency of the grackles that are temporarily held in aviaries in response to a novel environment and a novel object. Exploratory tendency is measured as the latency to approach and the closest approach distance to a novel object and a novel environment. Exploration assays occur twice for each bird: once near the beginning of their aviary time (“time 1”) and once again approximately 6 weeks later (“time 2”). Habituation may occur between time 1 and time 2, decreasing the novelty of the experimental setup. However, it is common practice to use the same setup across the repeated assays because it is very difficult to predict how threatening a novel object will be to a grackle. Therefore, if we accidentally introduce objects that are much more or much less threatening across the two time periods, this could obscure our ability to determine whether there are consistent individual differences with regard to these particular novel objects. We will analyze whether behavioral responses during assays are repeatable within individuals and whether exploration of a novel environment correlates with exploration of a novel object, indicating they are measures of the same inherent trait. If the two exploration measures are consistent within individuals and correlate with each other, we will choose as the exploration score the variable with the most data. If the two measures do not correlate, we will include both as independent variables.
Dependent variables
P1-P2
Home range size (square meters): an estimate calculated using the autocorrelated-Gaussian reference function kernel density estimate (AKDE), which is the only estimate of home range that accounts for autocorrelation due to the small time period between each of our GPS locations (Noonan et al. 2019). This estimate consists of the area enclosing the GPS location points for an individual grackle during its normal activities.
Autocorrelation of step length (meters): measured as the standard deviation of step lengths (the distance between two sequential GPS points)
Autocorrelation of turning angle (degrees): measured as the standard deviation of turning angles
Spatial location preference: measured as the repeatability of grackle occurrence in a given cell of a 5 x 5 m grid array across the landscape
One model will be run for each dependent variable
Independent variables
P1 and P1 alternatives 1-4
Exploration of novel environment: Latency to approach up to 20cm of a novel environment (that does not contain food) set inside a familiar environment (that contains maintenance diet away from the object) - OR - closest approach distance to the novel environment (choose the variable with the most data)
Exploration of novel object: Latency to approach up to 20cm of an object (novel or familiar, that does not contain food) in a familiar environment (that contains maintenance diet away from the object) - OR - closest approach distance to the object (choose the variable with the most data)
Sex: Male or female
History: the number of days the individual was temporarily held in the aviaries before data collection on space use began (0 indicates the grackle was only ever in the wild)
The number of known breeding sites (shade trees, date palms, marsh vegetation (K. Johnson and Peer 2001)) within the home range of each individual (data collected as part of Logan et al. 2019)
The number of human food source areas (dumpsters, cat food bowls, outdoor restaurant seating areas and parking lots) within the home range of each individual (data collected as part of Logan et al. 2019)
Scaled mass index (Peig and Green 2009) as a measure of energetic condition
Maximum group size observed across each individual’s focal follows (data collected as part of Logan et al. 2019)
Season (breeding or non-breeding)
P2
Site: Whether the grackle was from our study population located on the edge of the range (Northern California), the center of the original range (Central America), or the center of the current expanding edge (Arizona).
Sex: Male or female
History: the number of days the individual was temporarily held in the aviaries before data collection on space use began (0 indicates the grackle was only ever in the wild)
Population density (number of grackles per square meter in each study area: Arizona, California, Central America)
Season (breeding or non-breeding)
D. ANALYSIS PLAN
We do not plan to exclude any data and if there are missing data (e.g. if a bird had to be released before collecting their data at time 2) these birds will be excluded from analyses requiring data from times 1 and 2. Analyses will be conducted in R [current version 4.0.3; R Core Team (2017)] and Stan [version 2.18; Carpenter et al. (2017)].
When calculating our dependent variables we are interested in all movements by grackles, therefore we will not exclude any outlier relocations collected during “normal daily activities” (Calenge 2011). “Normal daily activities” indicate that grackles are not engaging in behaviors that would artificially skew their space use, for example mobbing a predator or the experimenter, or behavior before sunrise or after sunset when they are at the roost. Outside of these circumstances, we will include all data to detect space use movements.
From the GPS point locations collected on each individual, we will calculate home range size using the akde function in the R packages ctmm (Calabrese, Fleming, and Gurarie 2016) and sf (Pebesma 2018). We will use a Bayesian model (detailed below) to estimate the following parameters: mean and dispersion (variance) of step lengths and turning angles for each bird on each daily track (Pacheco-Cobos et al. 2019). We will determine whether these parameters governing movement are stable or variable within individuals across days. A small variance would indicate there is low variability (high repeatability) in the daily movement behaviors of the individual. Moreover, we will determine whether grackles show individual differences in consistent use of habitat by overlaying a grid array across the landscape. We will then create matrices describing the number of times a grackle was observed in each cell on each day. High autocorrelation among daily matrices indicates an individual that frequents the same spatial locations across days.
We will quantify region-specific grackle density by fitting a hierarchical model that accounts for imperfect detection. Specifically, we will use the model developed by Amundson, Royle, and Handel (2014) which integrates data on time of detection and distance estimates to account for the probability a bird is available to be detected (pa), and the probability it is detected given it is available (pd), respectively (Nichols, Thomas, and Conn 2009). All parameters (density, pa, and pd) will be modeled as a function of region. Because we expect that vocalization rates will be greater for males than females, we will also model pa as a function of sex. We will extract the estimate of expected point-level density for each region and use this estimate as the covariate in our movement models.
We will then model the relationship between independent variables describing bird-specific data on performance in the exploration tasks and other covariates (as outlined above), and bird-specific movement parameters (e.g. home range size, autocorrelation in step-size and turning angle, repeatability of spatial location preferences) as our dependent variables. We will use linear models and we will ensure assumptions of normality are met by checking that the residuals from our fitted models are normally distributed.
We will additionally test whether time spent in captivity might alter the social behavior of grackles when subsequently released into the wild by testing maximum group size observed across each individual’s focal follows as a function of the individual’s captivity history (the number of days the individual was temporarily held in the aviaries before data collection on space use began).
Ability to detect actual effects
To assess the effects that we will be able to detect given the expected sample size, we used G*Power [v.3.1; Faul et al. (2007); Faul et al. (2009)] to conduct power analyses based on confidence intervals. G*Power uses pre-set drop down menus and we chose the options that were as close to our analysis methods as possible (listed in each analysis below). These power analyses are not fully aligned with our study design, and the expected effect sizes are difficult to estimate due to the lack of prior data on this species; yet we are unaware of better options at this time.
Calculating home range size
#Load packages
library(ctmm)
library(sf)
#Point to the correct data file and load it
setwd("~/Documents/Grackle project/Space use")
pts<-read.csv("gtgr_points.csv", header = T)
#Set variables to ensure the models read the data properly
pts$Latitude = as.numeric(as.character(pts$Latitude))
pts$Longitude = as.numeric(as.character(pts$Longitude))
pts$Date = as.character(pts$Date)
pts$Time = as.character(pts$Time)
pts$da.ti = as.POSIXct(paste(pts$Date, pts$Time), format="%Y-%m-%d %H:%M:%S")
### Home range size can only be calculated when a bird has 20 or more relocations
### Count the number of relocations for each bird, then exclude individuals with less than 20
tmp = pts
tmp$count = 1
tmp = aggregate(count ~ Bird.Name, FUN = "sum", data = tmp)
pts.min = merge(pts, tmp, by = "Bird.Name", all = T)
pts.min = pts.min[-which(pts.min$count<20),]
pts.min$Bird.Name = as.factor(as.character(pts.min$Bird.Name))
#### Calculating home range size ####
#Turn dataframe into telemetry object
telem = pts.min[,c(1,2,4,5,7)] #only need bird name, date, lat and long, time
telem$timestamp = paste(telem$Date,telem$Time, sep=' ') #telem format requires date/time merged together for timestamp
telem = telem[,-c(2,5)] #remove other date and time columns
colnames(telem) = c("individual.local.identifier","location.lat","location.long","timestamp") #use required column names
telem = as.telemetry(telem, timeformat = "%Y-%m-%d %H:%M:%S",timezone = "")
for(i in 1:length(telem)){
GUESS = ctmm.guess(telem$i, interactive = FALSE)
m = ctmm.fit(telem$i, GUESS)
hr = akde(telem$i, m)
est[i] = summary(hr)$CI[2]
results = rbind(results,est[i])
}
### Normalize data for repeatability and regression models
hist(log(hr$area))
### Repeatability of home range in breeding and non-breeding seasons
hr_rpt = rpt(log(area) ~ (1|BirdID), grname = "BirdID", data = hr, datatype = "Gaussian", nboot = 500, npermut = 500)
summary(hr_rpt)
### Make data sheet with the experimental exploration measures combined with the home range measures
exp <- read.csv("Explore_combined.csv", header = T)
hr_exp = merge(exp, hr, by = "id", all = T)
write.csv(hr_exp, "space_use.csv")Code to create functions for analyzing movement behaviors
All scripts and code are available at https://github.com/ctross/grackleator.
### gracklenomial.R
gracklenomial <- function(GrackBins){
model_dat_2=list(
Trips=dim(GrackBins)[1],
Bins=dim(GrackBins)[2],
GrackBins=GrackBins)
model_code_2 <- '
data{
int Trips;
int Bins;
int GrackBins [Trips,Bins];
}
parameters{
simplex[Bins] A;
real<lower=0,upper=1> B;
}
model{
vector[Bins] C;
A ~ dirichlet(rep_vector(1,Bins));
B ~ beta(1,1);
for( i in 2:Trips){
C = to_vector(GrackBins[i-1])/sum(GrackBins[i-1]);
GrackBins[i] ~ multinomial( A*(1-B) + B*C );
}
}
'
m2 <- stan( model_code=model_code_2, data=model_dat_2,refresh=1,chains=1, control = list(adapt_delta = 0.9, max_treedepth = 13))
return(m2)
}
### gracklepars.R
gracklepars <- function(m1){
MuD_M<-data.frame(Group="Step-Size",Group2="Mean",Group3="Mean",Value=rstan::extract(m1,pars="MuD_M")$MuD_M) # Mean Step-Size over Days
DispD_M<-data.frame(Group="Step-Size",Group2="Mean",Group3="Dispersion",Value=rstan::extract(m1,pars="DispD_M")$DispD_M) # Dispersion in Mean Step-Size over days
MuD_D <-data.frame(Group="Step-Size",Group2="Dispersion",Group3="Mean",Value=rstan::extract(m1,pars="MuD_D")$MuD_D)# Mean Dispersion in Step-Size
DispD_D <-data.frame(Group="Step-Size",Group2="Dispersion",Group3="Dispersion",Value=rstan::extract(m1,pars="DispD_D")$DispD_D)# Dispersion in Dispersion in Step-Size
MuA_M <-data.frame(Group="Heading Change",Group2="Mean",Group3="Mean",Value=rstan::extract(m1,pars="MuA_M")$MuA_M)# Mean Angle Change over Days
DispA_M <-data.frame(Group="Heading Change",Group2="Mean",Group3="Dispersion",Value=rstan::extract(m1,pars="DispA_M")$DispA_M)# Dispersion in Mean Angle Change over days
MuA_D <-data.frame(Group="Heading Change",Group2="Dispersion",Group3="Mean",Value=rstan::extract(m1,pars="MuA_D")$MuA_D)# Mean Dispersion in Angle Change
DispA_D <-data.frame(Group="Heading Change",Group2="Dispersion",Group3="Dispersion",Value=rstan::extract(m1,pars="DispA_D")$DispA_D)# Dispersion in Dispersion in Angle Change
df_allx<-rbind(MuD_M,DispD_M,MuD_D,DispD_D,MuA_M,DispA_M,MuA_D,DispA_D)
ggplot(df_allx, aes(x=Value)) + geom_density(aes(group=Group3, colour=Group3, fill=Group3), alpha=0.3) + facet_wrap(Group2~Group,scales="free")+
labs(y="Regression parameters") + theme(strip.text.x = element_text(size=14,face="bold"),axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"))
}
### grackletrip.R
grackletrip <- function(m1){
m1a<-rstan::extract(m1,pars="AlphaAngle")
sample_eff<-apply(m1a$AlphaAngle,2,quantile,probs=c(0.05,0.5,0.95))
df_angle<-data.frame(Trip=c(1:dim(m1a$AlphaAngle)[2]),Group="Heading Change",Group2="Mean",
LI=sample_eff[1,],
Median=sample_eff[2,],
HI=sample_eff[3,])
m1d<-rstan::extract(m1,pars="AlphaDist")
sample_eff<-apply(m1d$AlphaDist,2,quantile,probs=c(0.05,0.5,0.95))
df_dist<-data.frame(Trip=c(1:dim(m1d$AlphaDist)[2]),Group="Step-Size",Group2="Mean",
LI=sample_eff[1,],
Median=sample_eff[2,],
HI=sample_eff[3,])
df_all1<-rbind(df_angle,df_dist)
m1a<-rstan::extract(m1,pars="DAngle")
sample_eff<-apply(exp(m1a$DAngle),2,quantile,probs=c(0.05,0.5,0.95))
df_angle<-data.frame(Trip=c(1:dim(m1a$DAngle)[2]),Group="Heading Change",Group2="Dispersion",
LI=sample_eff[1,],
Median=sample_eff[2,],
HI=sample_eff[3,])
m1d<-rstan::extract(m1,pars="SDDist")
sample_eff<-apply(exp(m1d$SDDist),2,quantile,probs=c(0.05,0.5,0.95))
df_dist<-data.frame(Trip=c(1:dim(m1d$SDDist)[2]),Group="Step-Size",Group2="Dispersion",
LI=sample_eff[1,],
Median=sample_eff[2,],
HI=sample_eff[3,])
df_all2<-rbind(df_angle,df_dist)
df_all<-rbind(df_all1,df_all2)
ggplot(df_all,aes(x=Trip,y=Median))+geom_point()+
geom_linerange(aes(ymin=LI,ymax=HI))+facet_wrap(Group2~Group,scales="free")+
labs(y="Regression parameters") + theme(strip.text.x = element_text(size=14,face="bold"),axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold")) +
scale_x_continuous(breaks= pretty_breaks())
}
### grackleate.R
grackleate <- function(AlphaDist=2.8,AlphaAngle=0,BetaDist=rep(0,10),BetaAngle=rep(0,10),SDDist=1,DAngle=2,Thresh=50, Lags=10, steps=150, seed=123456, Reps=30,N_Patches=13){
################################################################################# Set up Patches of Food
X_Mushrooms<-vector("list",N_Patches) # Location of prey
Y_Mushrooms<-vector("list",N_Patches) #
set.seed(seed) # Reset Seed
N_PerPatch<-rpois(N_Patches, 200) # Create some prey
for(i in 1:N_Patches){
X_Mushrooms[[i]]<-rnorm(N_PerPatch[i], runif(1,-2500,2500), rpois(1, 350))+100
Y_Mushrooms[[i]]<-rnorm(N_PerPatch[i], runif(1,-2500,2500), rpois(1, 350))+100
}
X_Mushrooms_per_step<-vector("list",steps) # Location of prey
Y_Mushrooms_per_step<-vector("list",steps) #
for(i in 1:steps){
X_Mushrooms_per_step[[i]]<-X_Mushrooms
Y_Mushrooms_per_step[[i]]<-Y_Mushrooms
}
################################################################################ Start Model
Store.X <- matrix(NA,ncol=Reps,nrow=steps-1)
Store.Y <- matrix(NA,ncol=Reps,nrow=steps-1)
for(q in 1:Reps){ # Run simulation several times
set.seed(seed) # Reset seed
N_PerPatch<-rpois(N_Patches, 200) # Choose number of items per patch
for(i in 1:N_Patches){
X_Mushrooms[[i]]<-rnorm(N_PerPatch[i], runif(1,-2500,2500), rpois(1, 350))+100 # Make some patches of prey
Y_Mushrooms[[i]]<-rnorm(N_PerPatch[i], runif(1,-2500,2500), rpois(1, 350))+100 #
}
X_Mushrooms_per_step<-vector("list",steps) # Location of prey
Y_Mushrooms_per_step<-vector("list",steps) #
for(i in 1:steps){
X_Mushrooms_per_step[[i]]<-X_Mushrooms
Y_Mushrooms_per_step[[i]]<-Y_Mushrooms
}
loc.x<-c() # Location of forager
loc.y<-c() #
speed<-c() # Speed or step size
heading<-c() # Absolute heading
d.heading<-c() # Heading change
Hits<-c() # Binary vector of encounters
loc.x[1:Lags]<-rep(0,Lags) # Intialize vectors
loc.y[1:Lags]<-rep(0,Lags) #
heading[1:Lags]<-rep(0,Lags) #
speed[1:Lags]<-rep(0,Lags) #
Hits[1:Lags]<-rep(0,Lags) #
plot(loc.x,loc.y,typ="l",ylim=c(-3500,3500),xlim=c(-3500,3500)) # Plot prey items
for(i in 1:N_Patches){ #
points(X_Mushrooms[[i]],Y_Mushrooms[[i]], col="red",pch=".") #
} #
set.seed(q*103)
################################################################################ Now model forager movement
for(s in (Lags+1): (steps-1)){
X_Mushrooms <- X_Mushrooms_per_step[[s]]
Y_Mushrooms <- Y_Mushrooms_per_step[[s]]
PredDist <- AlphaDist; # First calculate mean prediction conditional on encounters
for(k in 1:Lags){ #
PredDist <- PredDist + BetaDist[k]*ifelse(Hits[s-k]>0,1,0); #
} #
R<- exp(rnorm(1,PredDist,SDDist)) # Then simulate a step distance.
PredAngle <- AlphaAngle; # Again calculate mean prediction conditional on encounters
for(k in 1:Lags){ #
PredAngle <- PredAngle + BetaAngle[k]*ifelse(Hits[s-k]>0,1,0); #
} #
Theta<- rbeta(1,inv_logit(PredAngle)*DAngle, # And then simulate a directional change
(1-inv_logit(PredAngle))*DAngle)*180*ifelse(runif(1,0,1)>0.5,1,-1) # The 180 shifts from the unit to the maximum absolute distance
# The ifelse just chooses a random direction
heading[s]<-(heading[s-1]+Theta)%%360 # Store new heading. Note that the %% is the mod operation to wrap around if needed
d.heading[s] <- abs(Theta)/180 # Also store just the delta heading
speed[s] <- R # And the speed slash step size
ynew <- R * sin(deg2rad(heading[s])) # Now convert polar to Cartesian, to get the offset for a new x and y pair
xnew <- R * cos(deg2rad(heading[s])) #
loc.x[s]<-loc.x[s-1]+xnew # Make the new x and y pair
loc.y[s]<-loc.y[s-1]+ynew #
############################################################################### Now check for an encounter
Scrap2<-c()
for(i in 1:N_Patches){ # For each patch
Scrap<-rep(NA,length(X_Mushrooms[[i]]))
for(j in 1:length(X_Mushrooms[[i]])){ # For each mushroom in patch
Scrap[j]<- dist(loc.x[s],X_Mushrooms[[i]][j],loc.y[s],Y_Mushrooms[[i]][j]); # Calculate the distance from the forager to the mushroom
if(Scrap[j]<Thresh){ # If the forager is closer than the visual radius slash encounter threshold
points(X_Mushrooms[[i]][j],Y_Mushrooms[[i]][j], col="blue",pch=20) # Plot a hit
for(allgone in s: (steps-1)){
X_Mushrooms_per_step[[allgone]][[i]][j]<-99999 # If this run is for destructive foraging
X_Mushrooms_per_step[[allgone]][[i]][j]<-99999 # disappear food for all future timesteps
}
}}
Scrap2[i]<- ifelse(sum(ifelse(Scrap<Thresh,1,0))>0,sum(ifelse(Scrap<Thresh,1,0)),0) # Check for hits, also can replace sum(ifelse(Scrap<Thresh,1,0)) with 1
}
Hits[s]<-ifelse(sum(Scrap2)==0,0,sum(Scrap2)) # If there is a hit, then set encounters to 1
lines(loc.x,loc.y,ylim=c(-2500,2500),xlim=c(-2500,2500)) # Plot updates to the foragers path
}
Store.X[,q] <- loc.x
Store.Y[,q] <- loc.y
}
return(list(X=Store.X[(Lags+1):length(loc.x),],Y=Store.Y[(Lags+1):length(loc.y),]))
}
### grackleations.R
grackleations <- function(m2){
df2<-data.frame(Parameter="B",Value=rstan::extract(m2,pars="B")$B) # Mean Step-Size over Days
ggplot(df2, aes(x=Value)) + geom_density(aes(fill=Parameter), alpha=0.3) +
theme(strip.text.x = element_text(size=14,face="bold"),axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"))
}
### gracklebinner.R
gracklebinner <- function(tracks,nbin=c(15,15),ab_override=NULL){
Trips <- dim(tracks$X)[2]
GrackBins <- matrix(NA,nrow=Trips,ncol=nbin[1]*nbin[2])
bins<-bin2(cbind(c(z$X),c(z$Y)),nbin=nbin)
if(length(dim(ab_override))==0){
for( i in 1:Trips){
GrackBins[i,] <- c(bin2(cbind(z$X[,i],z$Y[,i]),nbin=nbin,ab=bins$ab)$nc)
}
} else{
for( i in 1:Trips){
GrackBins[i,] <- c(bin2(cbind(z$X[,i],z$Y[,i]),nbin=nbin,ab=ab_override)$nc)
}
}
return(GrackBins)
}
### grackleize.R
grackleize <- function(StoreX,StoreY) {
############################################################ Prepare Data
Trip <- dim(StoreX)[2]
MaxTicks<-dim(StoreX)[1]
Distance <-matrix(NA,ncol=max(Trip),nrow=MaxTicks)
Ang <-matrix(NA,ncol=max(Trip),nrow=MaxTicks)
AngDiff <-matrix(NA,ncol=max(Trip),nrow=MaxTicks)
N<-c()
for(j in 1:max(Trip)){
X <- StoreX[,j]
Y <- StoreY[,j]
N[j] <- length(X)
Dist<-c()
ang1<-c()
angdif<-c()
for(i in 2:N[j]){
Dist[i]<- dist(X[i],X[i-1],Y[i],Y[i-1]);
ang1[i]<- ang(X[i],X[i-1],Y[i],Y[i-1])
if(i>3){
angdif[i]<- ang.dif(ang1[i],ang1[i-1])
}}
Distance[1:N[j],j]<-Dist
Ang[1:N[j],j]<-ang1
AngDiff[1:N[j],j]<-angdif
}
for(i in 1:MaxTicks){
for(j in 1:max(Trip)){
Distance[i,j] <- ifelse(Distance[i,j]==0,
runif(1,0,min(Distance[which(Distance != 0)],na.rm=TRUE)),
Distance[i,j])
}}
Distance[is.na(Distance)] <- 999999
EmptyAngle <- fitdistr( ScaleBeta(c(na.omit(c(AngDiff)))),start=list(shape1=1,shape2=1),"beta")$estimate
for(i in 1:MaxTicks){
for(j in 1:max(Trip)){
AngDiff[i,j] <- ifelse(is.nan(AngDiff[i,j]), rbeta(1,EmptyAngle[1], EmptyAngle[2])*pi,AngDiff[i,j])
}}
AngDiff <- ScaleBeta(AngDiff)
AngDiff[is.na(AngDiff)] <- 999999
AngDiff <- AngDiff[c(-1,-2,-3),]
Distance <- Distance[c(-1,-2,-3),]
N <- N-3
MaxTicks<-MaxTicks-3
########################################################## Extract data for Stan
model_dat=list(
Distance=Distance,
AngDiff=AngDiff,
N=N,
MaxTrip=max(Trip),
MaxTicks=MaxTicks)
################################################################# Fit Stan Model
model_code_1 <- '
data{
int MaxTrip;
int MaxTicks;
int N[MaxTrip];
real Distance[MaxTicks,MaxTrip];
real AngDiff[MaxTicks,MaxTrip];
}
parameters {
real MuD_M;
real MuA_M;
real<lower=0> DispD_M;
real<lower=0> DispA_M;
real MuD_D;
real MuA_D;
real<lower=0> DispD_D;
real<lower=0> DispA_D;
vector[MaxTrip] AlphaDist_Raw;
vector[MaxTrip] SDDist_Raw;
vector[MaxTrip] AlphaAngle_Raw;
vector[MaxTrip] DAngle_Raw;
}
transformed parameters{
vector[MaxTrip] AlphaDist;
vector[MaxTrip] SDDist;
vector[MaxTrip] AlphaAngle;
vector[MaxTrip] DAngle;
AlphaDist = MuD_M + DispD_M*AlphaDist_Raw;
SDDist = MuD_D + DispD_D*SDDist_Raw;
AlphaAngle = MuA_M + DispA_M*AlphaAngle_Raw;
DAngle = MuA_D + DispA_D*DAngle_Raw;
}
model{
MuD_M~normal(0,5);
MuA_M~normal(0,5);
DispD_M~cauchy(0,1);
DispA_M~cauchy(0,1);
MuD_D~normal(0,5);
MuA_D~normal(0,5);
DispD_D~cauchy(0,1);
DispA_D~cauchy(0,1);
AlphaDist_Raw~normal(0,1);
SDDist_Raw~normal(0,1);
AlphaAngle_Raw~normal(0,1);
DAngle_Raw~normal(0,1);
for(j in 1:MaxTrip){
{
vector[N[j]] PredAngle;
vector[N[j]] PredDist;
vector[N[j]] Dist;
vector[N[j]] AngleDifference;
for(i in 1:N[j]){
PredDist[i] = AlphaDist[j];
}
for(i in 1:N[j]){
PredAngle[i] = AlphaAngle[j];
}
for(i in 1:N[j]){
Dist[i] = Distance[i,j];
AngleDifference[i] = AngDiff[i,j];
}
Dist ~ lognormal(PredDist,exp(SDDist[j]));
AngleDifference ~ beta(inv_logit(PredAngle)*exp(DAngle[j]), (1-inv_logit(PredAngle))*exp(DAngle[j]));
}}
}
'
m1 <- stan( model_code=model_code_1, data=model_dat,refresh=1,chains=1, control = list(adapt_delta = 0.9, max_treedepth = 15))
return(m1)
}
### setup.R
# Load libraries
library(MASS)
library(mvtnorm)
library(fields)
library(ggplot2)
library(rethinking)
library(sfsmisc)
library(ash)
library(reshape)
library(maptools)
library(gridExtra)
library(scales)
library(msm)
library(maps)
library(grid)
library(xtable)
# Define functions
dist2<-function(a,b){return(sqrt( (b[2]-a[2])^2 + (b[1]-a[1])^2 ))}
dist<-function(x2a,x1a,y2a,y1a){return(sqrt((x2a-x1a)^2 + (y2a-y1a)^2))}
rad2deg <- function(rad) {(rad * 180) / (pi)}
deg2rad <- function(deg) {(deg * pi) / (180)}
ang.dif <- function(x,y) {min((2 * pi) - abs(x - y), abs(x - y))}
lp_dist <-function(a,b,c){
# Return distance between a point and line segment
t <- -(((a[1]-c[1])*(b[1]-a[1]) + (a[2]-c[2])*(b[2]-a[2])) / ((b[1]-a[1])^2 + (b[2] - a[2])^2 ))
if(t >0 & t <1){
numer <- abs( (b[2]-a[2])*c[1] -(b[1]-a[1])*c[2] + b[1]*a[2] - b[2]*a[1])
denom <- sqrt( (b[2]-a[2])^2 + (b[1]-a[1])^2 )
return(numer/denom)
}
else{
d1 <- dist2(a,c)
d2 <- dist2(b,c)
return( min(d1,d2))
}
}
ScaleBeta <- function(X){
a<-0
b<-pi
y<-X
Samp<-50
y2 <- (y-a)/(b-a)
y3<-(y2*(Samp - 1) + 0.5)/Samp
y3}
ang<-function(x2,x1,y2,y1){
Theta<-seq(0,2*pi,length.out=100)
DB<-cbind(cos(Theta),sin(Theta))
r <- dist(x2,x1,y2,y1)
x0<-(x2-x1)/r
y0<-(y2-y1)/r
(atan2(y0,x0) + pi )
}
################################################################# Set parameters
Thresh <- 50 # Threshold of visual range
Lags<-10 # Lags at which there are effects of encounters on movement, kinda hardcoded---careful with changes
Steps<-150 # Length of simulation
Reps <- 30 # Number of replications
# Adaptive model parameters
Phi0A <- 3.2
Psi0A <- 0
PhiA <- c(-0.76, -0.65, -0.58, -0.43, -0.29, -0.15, -0.03,0,0,0)
PsiA <- c(0.5, 0.3, 0.25, 0.2, 0.17, 0.15, 0.13, 0.11, 0.09, 0.08)
OmegaA <- 0.4
EtaA <- 2
# Levy parameters
Phi0L <- 2.8
Psi0L <- 0
PhiL <- c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
PsiL <- c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
OmegaL <- 1
EtaL <- 2
# Brownian parameters
Phi0B <- 25
Psi0B <- 0
PhiB <- c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
PsiB <- c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
OmegaB <- 15
EtaB <- 2Modeling bird movement behaviors: step length, turning angle, spatial location preference
#This is an R package for modeling bird movement (available at https://github.com/ctross/grackleator)
#Install by running on R:
library(devtools)
install_github("Ctross/grackleator")
#Some quick examples:
#First, lets look at movement dynamics:
# Load library and attach data
library(grackleator)
# Simulate tracks from 1 grackle over 30 trips
z <- grackleate(AlphaDist=3.8,AlphaAngle=0,SDDist=1.5,DAngle=2)
# Analyze results of trips for step-size and angle change
m1 <- grackleize(z$X,z$Y)
# Plot results over trips
grackletrip(m1)
# Plot bird-specific parameters
gracklepars(m1)
#And now habitat selection:
# Load library and attach data
library(grackleator)
# Simulate tracks from 1 grackle over 30 trips
z <- grackleate(AlphaDist=3.8,AlphaAngle=0,SDDist=1.5,DAngle=2)
# Bin data on 2D grid
GrackBins<- gracklebinner(z,ab_override=matrix(c(-3000,-3000,3000,3000),nrow=2,ncol=2))
# Model location data
m2 <- gracklenomial(GrackBins)
# Plot environmental hotspots
image(matrix(colSums(GrackBins), nrow=15,ncol=15)) # Overall
image(matrix(GrackBins[1,], nrow=15,ncol=15)) # Day 1
image(matrix(GrackBins[2,], nrow=15,ncol=15)) # Day 2
image(matrix(GrackBins[3,], nrow=15,ncol=15)) # Day 3
image(matrix(get_posterior_mean(m2,pars="A"), nrow=15,ncol=15)) # Overall
# Plot bird-specific parameters
grackleations(m2)
# Now create data with temporal correlations
GrackBins2 <- GrackBins
A <- get_posterior_mean(m2,pars="A")
A <- A/sum(A)
B <- 0.8
for(i in 2:30)
GrackBins2[i,] <- rmultinom(1,sum(GrackBins2[i-1,]),A*(1-B) + B*(GrackBins2[i-1,]/sum(GrackBins2[i-1,])))
# Analyze the new data
m3 <- gracklenomial(GrackBins2)
# Plot environmental hotspots
image(matrix(colSums(GrackBins2), nrow=15,ncol=15)) # Overall
image(matrix(GrackBins2[1,], nrow=15,ncol=15)) # Day 1
image(matrix(GrackBins2[2,], nrow=15,ncol=15)) # Day 2
image(matrix(GrackBins2[3,], nrow=15,ncol=15)) # Day 3
image(matrix(get_posterior_mean(m3,pars="A"), nrow=15,ncol=15)) # Post
# Plot bird-specific parameters
grackleations(m3)H1: P1 - Exploration measured in captivity relates to space use behavior
To roughly estimate our ability to detect actual effects, we ran a power analysis in G*Power with the following settings: test family=F tests, statistical test=linear multiple regression: Fixed model (R^2 deviation from zero), type of power analysis=a priori, alpha error probability=0.05. We reduced the power to 0.70 and increased the effect size until the total sample size in the output matched our projected sample size (n=60). The protocol of the power analysis is here:
Input:
Effect size f² = 0,24
α err prob = 0,05
Power (1-β err prob) = 0,7
Number of predictors = 8
Output:
Noncentrality parameter λ = 14,4000000
Critical F = 2,1260234
Numerator df = 8
Denominator df = 51
Total sample size = 60
Actual power = 0,7039242
This means that, with our minimum sample size of 60, we have a 70% chance of detecting a small effect (approximated at f2=0.20 by Cohen 1988).
data <- read.csv("Space_use.csv", header = T)
# Home range
#m1 = lm(log(area) ~ ExpObj + ExpEnv + Sex + History + Breeding + Feeding + Condition + Group + Season, data = data)
hist(m1$resid)
summary(m1)
# Step length
#m2 = lm(log(std_step)) ~ ExpObj + ExpEnv + Sex + History + Breeding + Feeding + Condition + Group + Season, data = data)
hist(m2$resid)
summary(m2)
# Turning angle
#m3 = lm(log(std_angle)) ~ ExpObj + ExpEnv + Sex + History + Breeding + Feeding + Condition + Group + Season, data = data)
hist(H3$resid)
summary(m3)
# Spatial preferences
#m4 = lm(log(loc_pref)) ~ ExpObj + ExpEnv + Sex + History + Breeding + Feeding + Condition + Group + Season, data = data)
hist(H4$resid)
summary(m4)H2: P2 - Space use behaviors vary among populations across the range
To roughly estimate our ability to detect actual effects, we ran a power analysis in G*Power in the same way as for Hypothesis 1. The protocol of the power analysis is here:
Input:
Effect size f² = 0,18
α err prob = 0,05
Power (1-β err prob) = 0,7
Number of predictors = 4
Output:
Noncentrality parameter λ = 10,6200000
Critical F = 2,5429175
Numerator df = 4
Denominator df = 54
Total sample size = 59
Actual power = 0,7025819
This means that, with our minimum sample size of 60, we have a 70% chance of detecting a small effect (approximated at f2=0.20 by Cohen 1988).
data <- read.csv("Space_use.csv", header = T)
# Home range
#m1 = lm(log(area) ~ Site + Sex + History + Season, data = data)
hist(m1$resid)
summary(m1)
# Step length
#m2 = lm(log(std_step)) ~ Site + Sex + History + Season, data = data)
hist(m2$resid)
summary(m2)
# Turning angle
#m3 = lm(log(std_angle)) ~ Site + Sex + History + Season, data = data)
hist(m3$resid)
summary(m3)
# Spatial preference
#m4 = lm(log(loc_pref)) ~ Site + Sex + History + Season, data = data)
hist(m4$resid)
summary(m4)E. ETHICS
This research is carried out in accordance with permits from the:
- US Fish and Wildlife Service (scientific collecting permit number MB76700A-0,1,2)
- US Geological Survey Bird Banding Laboratory (federal bird banding permit number 23872)
- Arizona Game and Fish Department (scientific collecting license number SP594338 [2017], SP606267 [2018], SP639866 [2019], and SP402153 [2020])
- Institutional Animal Care and Use Committee at Arizona State University (protocol number 17-1594R)
G. FUNDING
This research is funded by the Department of Human Behavior, Ecology and Culture at the Max Planck Institute for Evolutionary Anthropology.
I. CONFLICT OF INTEREST DISCLOSURE
We, the authors, declare that we have no financial conflicts of interest with the content of this article. Corina Logan is a Recommender and on the Managing Board at PCI Ecology.
J. ACKNOWLEDGEMENTS
We thank Melissa Wilson Sayres for sponsoring our affiliations at Arizona State University (ASU); Kevin Langergraber for serving as the local PI on the ASU IACUC; Kristine Johnson for technical advice on great-tailed grackles; Julia Cissewski and Sophie Kaube for tirelessly solving problems involving financial transactions and contracts; Richard McElreath for project support; and our research assistants: Aldora Messinger, Elysia Mamola, Michael Guillen, Rita Barakat, Adriana Boderash, Olateju Ojekunle, August Sevchik, Justin Huynh, Jennifer Berens, Michael Pickett, Amanda Overholt, Emily Blackwell, Kaylee Delcid, Brynna Hood, Sam Bowser, Sierra Planck, Samuel Munoz, Sawyer Lung.