Load Data and Wrangle
# This function makes sure all circular values are coded as circular, and run as -pi to pi
wranglecircular=function(Angles){
Angles = as.circular(Angles, type="angles",units="radians",
template="none",modulo="asis",
zero=0,rotation="counter")
Angles=(Angles+(2*pi))%%(2*pi)
Angles=if_else(Angles<pi,Angles,Angles-(2*pi))
Angles
}
# Load The Trial Summary Data
## Colmuns:
## GS = The group size
## Trial = The trial number (the NA for GS 1 is Trial 14)
## Center_X/Y = The x,y coordinates of the center of the pan
## Beach_X/Y = The x,y coordinates of the spot on the edge of the pan closest to the beach/wrackline
## AngleBeach = The angle between the center of the pan and the beach
## Date = Date of the trial in Year-Month-Day.
## Weather = A qualitative descriptor of the weather
## Temp = Temperature at the start of the trial in Fahrenheit
## hour/minutes = The hour and minute of the day when the trial began in 24hr
## Dig = The name of the person who digitized the trial
## Dist_Actual = The sum of the distance traveled between each step (cm)
## Dist_StrtLine = The beeline distance between their starting and ending coordinates (cm)
## Sinuosity = Dist_StrtLine/Dist/Actual. Not an accurate measurement of path wiggliness b/c of non-standard steplengths
## Time = The total number of seconds between when the animal started moving to when they finished the trial
## FinalPosition = The angle between where the animal ended the trial and the center of the arena relative to the wrackline. Renamed to Final Orientation
## HeadingMean = The circular mean of each heading, weighted by the distance traveled relative to the wrackline
## StartHeading = The bearing of their first step, relative to the wrackline. Renamed to InitialBearing
## Velocity = Dist_Actual / Time in cm per s
## Start_X/Y = Their starting X,Y location in the original coordinate system
## Start_X/Y.T = Their starting X,Y location transformed to a new coordinate system where the center of the pan is at 0,0, and the angle between the center of the pan and the beach is = 0 (rotated such that the beach is always to the right of all figures)
## StartLocation = Their angle between where the animal started the trial and the center of the arena relative to the wrackline. Not equivalent to StartHeading which is the angle between their first and second step. Renamed to InitialOrientation
## End_X/Y = Their final X,Y location in the original coordinate system
## Latency = The time between when the trial started (the cup was lifted) and when the animal began moving in s. This time is not included in their total time
## Animal = a designator for which animal was being digitized within each trial
## Dist_Diff = Dist_Actual-Dist_StrtLine. Renamed to ExcessDist
## TOD = Time of Day at the beginning of the trial
## Smoky = A designator for the days during which there was a nearby forest fire
## Year = The year the trial was conducted
## Latency2 = log(Latency), used for fig S5
## GS.Trial = GS and Trial combined. Used for MEM's
## FinalPosition_pan = The angle between where the animal ended the trial and the center of the arena relative to the arena. Renamed to FinalOrientation_pan
## FinalPosition/StartHeading.L = A linearized version of the Final position and Start Heading. Renamed to InitialBearing.L and FinalOrientation.L
GSall.alltrials.Outlier=read.csv("GSall.csv") |>
mutate(
AngleBeach=wranglecircular(AngleBeach),
FinalOrientation=wranglecircular(FinalPosition),
HeadingMean=wranglecircular(HeadingMean),
InitialBearing=wranglecircular(StartHeading),
InitialOrientation=wranglecircular(StartLocation),
FinalOrientation.Pan=wranglecircular(FinalPosition.Pan)
) |>
rename(InitialBearing.L=StartHeading.L,
FinalOrientation.L=FinalPosition.L,
ExcessDist=Dist_Diff)
# Remove outliers
GSall.alltrials=GSall.alltrials.Outlier |>
filter(Time<18)
# A df where MRV is calculated for each trial (unable to be calculated for GS1 as there was only one animal). Used for MEM's.
GSall.alltrials.mrv=GSall.alltrials |>
group_by(GS, Trial, GS.Trial, Weather, Temp, TOD, Smoky, Year) |>
summarise(FinalMRV = rho.circular(FinalOrientation), InitialMRV = rho.circular(InitialBearing)) |>
filter(GS != 1)
# Load the Full Track Data
## Columns:
## X/Y = x,y of each step in the original coordinate system from the moment the animal starts moving till the trial ends
## Center_X/Y (changed to X/Y.Centered) = x,y of each step rel. to CENTER OF ARENA (in this coordinate system, center=0,0 but the beach is still AngleBeach)
## Dist_Center = distance between each step to CENTER OF ARENA (calculated after centering each x,y)
## AngleCenter = angle between each step and CENTER OF ARENA (calculated after centering each x,y)
## MinPrev_X/Y = x,y of each step rel. to PREVIOUS STEP (calculated on the original x,y)
## Dist_Each = distance between EACH STEP (calculated on the original x,y)
## AngleCenter_Transformed = angle between each step and CENTER OF ARENA, transformed such that 0=beach (calculated as AngleCenter - AngleBeach)
## X/Y.Transformed = new x,y locations transfored (in this coordinate system, center=0,0 and the angle towards the beach is = 0. Calculated as Dist_Center*cos or sin of AngleCenter_Transformed)
## Heading = the individuals actual bearing, relative to the beach. 0 = towards the beach (calculated using X/Y.Transformed)
## Time = frame of that particular step
GSall.tracks=read.csv("GSalltracks.csv") |>
mutate(
Heading=wranglecircular(Heading),
AngleCenter_Transformed=wranglecircular(AngleCenter_Transformed)
)
# Load the Bootstrapped Data (Bootstrap occurs in markdown document "Nonparametric Bootstrap Analysis.Rmd" but because it takes some time to run we saved the outputs from that file and load them here.)
GS1.5=read.csv("Bootstrapped Underlying Dist for GS5.csv") |>
mutate(
InitialBearing=wranglecircular(StartHeading),
FinalOrientation=wranglecircular(FinalPosition),
InitialOrientation=wranglecircular(StartPosition)
) |>
rename(InitialMRV="rho_StartHeading",
FinalMRV="rho_FinalPosition",
ExcessDist="Dist_Diff")
GS1.10=read.csv("Bootstrapped Underlying Dist for GS10.csv") |>
mutate(
InitialBearing=wranglecircular(StartHeading),
FinalOrientation=wranglecircular(FinalPosition),
InitialOrientation=wranglecircular(StartPosition)
) |>
rename(InitialMRV="rho_StartHeading",
FinalMRV="rho_FinalPosition",
ExcessDist="Dist_Diff")
df.GS1.mean=read.csv("Averaged GS1 Values.csv") |>
mutate(
InitialBearing=wranglecircular(StartHeading),
FinalOrientation=wranglecircular(FinalPosition),
InitialOrientation=wranglecircular(StartPosition)
) |>
rename(InitialMRV="rho_StartHeading",
FinalMRV="rho_FinalPosition",
ExcessDist="Dist_Diff")
Main Text Figures and Calculations
Mean angle and Confidence Interval for Final Orientation Relative to
the Beach
GSall.alltrials |>
#Filter to the data we want
filter(Weather=="Sunny", Smoky=="Not Smoky", Time<18) |>
#Looking at Final Orientation
rename(Variable=FinalOrientation) |>
#calculate the mean per trial
group_by(GS.Trial) |>
summarise(Variable=mean.circular(Variable)) |>
#convert to radians, -pi to pi, where 0 = towards the beach
mutate(
Variable=(Variable+(2*pi))%%(2*pi),
Variable=if_else(Variable<pi,Variable,Variable-(2*pi))
) |>
ungroup() |>
#calculate values over all trials
summarise(MeanHeading=mean.circular(Variable), #mean angle
Rho=rho.circular(Variable), #rho
pvalue=round(rayleigh.test(Variable)$p.value,3), #pvalue
CI1=mle.vonmises.bootstrap.ci(Variable)$mu.ci[1], #95% confidence interval on the mean angle (upper angle)
CI2=mle.vonmises.bootstrap.ci(Variable)$mu.ci[2], #(lower angle)
FO.CI=mean.circular(abs(MeanHeading-CI1),abs(MeanHeading-CI2) ) #the difference between the mean angle and the confidence interval angles
) |>
#convert back to degrees
mutate(MeanHeading=MeanHeading*(180/pi),
FO.CI=FO.CI*(180/pi))
## # A tibble: 1 × 6
## MeanHeading Rho pvalue CI1 CI2 FO.CI
## <circular> <dbl> <dbl> <circular> <circular> <circular>
## 1 17.86063 0.817 0 0.1888232 0.4198099 7.041857
Table 2: The mean and sd for each movement variable for each group
size
(Table2= GSall.alltrials |>
#filter to the data we want
filter(Time<18, Weather=="Sunny", Smoky=="Not Smoky") |>
group_by(GS) |>
#rename Group Sizes and set as factor
mutate(GS=recode(GS,"1"="GS1","5"="GS5", "10"="GS10"))|>
mutate(GS=factor(GS,levels=c("GS1","GS5","GS10"))) |>
#calculate the mean and standard deviation of all variables of interest by group size
summarise(
Latency.mean=mean(Latency),Latency.sd=sd(Latency), #latency
Time.mean=mean(Time),Time.sd=sd(Time), #Total Trial Time
Speed.mean=mean(Velocity),Speed.sd=sd(Velocity), #Speed
ExcessDist.mean=mean(ExcessDist),ExcessDist.sd=sd(ExcessDist), #Excess Distance
InitialBearing.mean=mean.circular(InitialBearing),
IB.sd=sd.circular(InitialBearing), #Inital Bearing
FinalOrientation.mean=mean.circular(FinalOrientation),
FO.sd=sd.circular(FinalOrientation) #Final Orientation
) |>
#convert angles into degrees
mutate(
InitialBearing.mean=(InitialBearing.mean)*(180/pi),
FinalOrientation.mean=FinalOrientation.mean*(180/pi),
IB.sd=(IB.sd)*(180/pi),
FO.sd=(FO.sd)*(180/pi),
) |>
#Get the MRV for the initial bearing and final orientation from bootstrapped data
mutate(InitialMRV=c(df.GS1.mean$InitialMRV,mean(GS1.5$InitialMRV),mean(GS1.10$InitialMRV))) |>
mutate(FinalMRV=c(df.GS1.mean$FinalMRV,mean(GS1.5$FinalMRV),mean(GS1.10$FinalMRV)))|>
#round values to make the table legible
mutate(across(2:13, \(x) round(x, 1))) |>
mutate(across(14:15, \(x) round(x, 2))) |>
ungroup() |>
# Add standard deviation to mean calculation with a +/-
reframe(GS=GS,
"Latency (s)"=paste0(Latency.mean," ± ",Latency.sd),
"Time (s)"=paste0(Time.mean," ± ",Time.sd),
"Speed (cm/s)"=paste0(Speed.mean," ± ",Speed.sd),
"Excess Dist (cm)"=paste0(ExcessDist.mean," ± ",ExcessDist.sd),
"Initial Bearing (deg.)"=paste0(InitialBearing.mean," ± ",IB.sd),
"Initial MRV"=as.character(InitialMRV),
"Final Orientation (deg.)"=paste0(FinalOrientation.mean," ± ",FO.sd),
"Final MRV"=as.character(FinalMRV)
) |>
#pivot to correct orientation
pivot_longer(cols=-GS, names_to="Metrics") |>
pivot_wider(names_from="GS", values_from="value")
)
## # A tibble: 8 × 4
## Metrics GS1 GS5 GS10
## <chr> <chr> <chr> <chr>
## 1 Latency (s) 12.4 ± 6.6 1.3 ± 2.8 1.2 ± 3.4
## 2 Time (s) 3.8 ± 2.2 3.8 ± 1.7 4.1 ± 2.4
## 3 Speed (cm/s) 5.9 ± 1.7 5.3 ± 2.3 4.9 ± 2.4
## 4 Excess Dist (cm) 8.5 ± 7.4 7.1 ± 5.4 5.8 ± 4.7
## 5 Initial Bearing (deg.) -54.9 ± 121.6 5.7 ± 86.8 16.6 ± 82.7
## 6 Initial MRV 0.11 0.34 0.37
## 7 Final Orientation (deg.) 19.2 ± 54.5 16.7 ± 40.3 18.6 ± 44.6
## 8 Final MRV 0.64 0.78 0.74
# Using the code from above but without GS5 or the SD's, to calculate the difference between the GS1 and GS10 values for numbers quoted in text
GSall.alltrials |>
filter(Time<18, Weather=="Sunny", Smoky=="Not Smoky",GS!=5) |>
group_by(GS) |>
mutate(GS=recode(GS,"1"="GS1", "10"="GS10"))|>
mutate(GS=factor(GS,levels=c("GS1","GS10"))) |>
summarise(
Latency.mean=mean(Latency),
Time.mean=mean(Time),
Speed.mean=mean(Velocity),
ExcessDist.mean=mean(ExcessDist),
InitialBearing.mean=mean.circular(InitialBearing),
FinalOrientation.mean=mean.circular(FinalOrientation)
) |>
mutate(
InitialBearing.mean=as.numeric((InitialBearing.mean)*(180/pi)),
FinalOrientation.mean=as.numeric(FinalOrientation.mean*(180/pi) )) |>
mutate(InitialMRV=c(df.GS1.mean$InitialMRV,mean(GS1.10$InitialMRV))) |>
mutate(FinalMRV=c(df.GS1.mean$FinalMRV,mean(GS1.10$FinalMRV)))|>
mutate(across(2:7, \(x) round(x, 1))) |>
mutate(across(8:9, \(x) round(x, 2))) |>
ungroup() |>
pivot_longer(cols=-GS, names_to="Metrics") |>
pivot_wider(names_from="GS", values_from="value")|>
mutate(diff10=GS10-GS1)
## # A tibble: 8 × 4
## Metrics GS1 GS10 diff10
## <chr> <dbl> <dbl> <dbl>
## 1 Latency.mean 12.4 1.2 -11.2
## 2 Time.mean 3.8 4.1 0.300
## 3 Speed.mean 5.9 4.9 -1
## 4 ExcessDist.mean 8.5 5.8 -2.7
## 5 InitialBearing.mean -54.9 16.6 71.5
## 6 FinalOrientation.mean 19.2 18.6 -0.600
## 7 InitialMRV 0.11 0.37 0.26
## 8 FinalMRV 0.64 0.74 0.1
Supplemental Information Figures and Analyses
Supplemental Table S1: Mixed Effect Model results for Cloud Cover
and Year v.s. the eight movement variables of interest
#load function to make the MEM's
source("Mixed Effect Model Function.R")
# Cloud Cover
CloudCover=new.lmer("Latency","Weather",randomeffects = 2)
CloudCover=rbind(CloudCover, new.lmer("Time","Weather",randomeffects = 2))
CloudCover=rbind(CloudCover, new.lmer("Velocity","Weather",randomeffects = 2))
CloudCover=rbind(CloudCover, new.lmer("ExcessDist","Weather",randomeffects = 2))
CloudCover=rbind(CloudCover, new.lmer("InitialBearing.L","Weather",randomeffects = 2))
CloudCover=rbind(CloudCover, new.lmer("InitialMRV","Weather",randomeffects = 1,df = GSall.alltrials.mrv))
CloudCover=rbind(CloudCover, new.lmer("FinalOrientation.L","Weather",randomeffects = 2))
CloudCover=rbind(CloudCover, new.lmer("FinalMRV","Weather",randomeffects = 1,df = GSall.alltrials.mrv))
#Year
Year=new.lmer("Latency","Year",randomeffects = 3)
Year=rbind(Year, new.lmer("Time","Year",randomeffects = 3))
Year=rbind(Year, new.lmer("Velocity","Year",randomeffects = 3))
Year=rbind(Year, new.lmer("ExcessDist","Year",randomeffects = 3))
Year=rbind(Year, new.lmer("InitialBearing.L","Year",randomeffects = 3))
Year=rbind(Year, new.lmer("InitialMRV","Year",randomeffects = 4, GSall.alltrials.mrv))
Year=rbind(Year, new.lmer("FinalOrientation.L","Year",randomeffects = 3))
Year=rbind(Year, new.lmer("FinalMRV","Year",randomeffects = 4, GSall.alltrials.mrv))
#Bind together
(TableS1=rbind(CloudCover, Year) |>
relocate(Main.Effect,.before = Response) |>
mutate(Estimate=round(Estimate,5),
Std.Error=round(Std.Error,5)) |>
mutate(Response=recode(Response,"Velocity"="Speed"))
)
## Main.Effect Response Estimate Std.Error pvalue SigCode
## 1 Weather Latency 0.44909 0.34581 0.192 n.s.
## 2 Weather Time -1.15480 0.28191 0.000 ***
## 3 Weather Speed 0.34897 0.24408 0.136 n.s.
## 4 Weather ExcessDist -3.00548 0.80733 0.000 ***
## 5 Weather InitialBearing.L -0.12843 0.09805 0.189 n.s.
## 6 Weather InitialMRV 0.06635 0.05189 0.181 n.s.
## 7 Weather FinalOrientation.L -0.61067 0.08277 0.000 ***
## 8 Weather FinalMRV 0.21877 0.03807 0.000 ***
## 9 Year Latency 0.26989 0.18165 0.155 n.s.
## 10 Year Time 0.20309 0.17095 0.313 n.s.
## 11 Year Speed -0.19949 0.14659 0.267 n.s.
## 12 Year ExcessDist -0.61855 0.47708 0.154 n.s.
## 13 Year InitialBearing.L -0.09129 0.05108 0.073 n.s.
## 14 Year InitialMRV -0.00509 0.02897 0.980 n.s.
## 15 Year FinalOrientation.L -0.11226 0.05043 0.027 n.s.
## 16 Year FinalMRV 0.03996 0.02267 0.076 n.s.
Supplemental Figure S7: The histogram plots from the bootstrap
analysis withh associated p-values for all movement variables for GS5
and GS10
#load functions to make the plots for the linear and circular variables
source("Bootstrap Histogram Plot Function.R")
source("Circular Bootstrap Histogram Plot Function.R")
# make the plots for each of the variables of interest, comparing GS1 to both GS5 and GS10
L5=bootstrap.histogram.plot("Latency","GS5", bins=70,maxy = 0.15, title="", xaxis = "Latency (s)")+scale_x_continuous(breaks = c(0,5,10))#+theme(plot.title = element_text(hjust=1.5))
L10=bootstrap.histogram.plot("Latency","GS10", bins=60, maxy = 0.15,yaxis = F, title = "",xaxis = "Latency (s)")+scale_x_continuous(breaks = c(0,5,10))
V5=bootstrap.histogram.plot("Velocity","GS5",maxy=0.15,yaxis = T, title="", xaxis = "Speed (cm/s)")
V10=bootstrap.histogram.plot("Velocity","GS10",maxy=0.15,yaxis = F, title="", xaxis = "Speed (cm/s)")
DD5=bootstrap.histogram.plot("ExcessDist","GS5", maxy=0.15, yaxis=F,title = "",xaxis = "Excess Dist (cm)")
DD10=bootstrap.histogram.plot("ExcessDist","GS10",maxy=0.15,yaxis = F, title="",xaxis = "Excess Dist G(cm)")
T5=bootstrap.histogram.plot("Time","GS5", maxy=0.15, yaxis = F, title = "", xaxis = "Time (s)")
T10=bootstrap.histogram.plot("Time","GS10", maxy=0.15, yaxis = F, title="",xaxis = "Time (s)")
RF5=bootstrap.histogram.plot("FinalMRV","GS5",maxy=0.15,yaxis = T,title = "", xaxis = "Final MRV")
RF10=bootstrap.histogram.plot("FinalMRV","GS10",maxy=0.15,yaxis = F, title = "", xaxis = "Final MRV")
RS5=bootstrap.histogram.plot("InitialMRV","GS5",maxy=0.15,yaxis = F,title = "", xaxis="Initial MRV")
RS10=bootstrap.histogram.plot("InitialMRV","GS10",maxy=0.15,yaxis = F, title="", xaxis="Initial MRV")
SH5.C=bootstrap.histogram.plot.circular("InitialBearing", GS1.5,"#AF4E4E",title = "",xaxis="Initial Bearing",R=500)
SH10.C=bootstrap.histogram.plot.circular("InitialBearing", GS1.10, "#FBAD4F",title = "",xaxis="Initial Bearing",R=500)
F5.C=bootstrap.histogram.plot.circular("FinalOrientation", GS1.5,"#AF4E4E",R=500,title = "",xaxis="Final Orientation")
F10.C=bootstrap.histogram.plot.circular("FinalOrientation", GS1.10, "#FBAD4F",R=500,title = "",xaxis="Final Orientation")
p5=L5+T5+V5+DD5+RF5+RS5+F5.C+SH5.C+plot_layout(ncol=2,nrow=4, heights=c(5,5,5,5), widths=c(4,4))+plot_annotation(title="GS1 vs. GS5", theme=theme(plot.title = element_text(size = 15,hjust=0.5)))
p10=L10+T10+V10+DD10+RF10+RS10+F10.C+SH10.C+plot_layout(ncol=2,nrow=4, heights=c(5,5,5,5), widths=c(4,4))+plot_annotation(title="GS1 vs. GS10", theme=theme(plot.title = element_text(size=15,hjust=0.5)))
wrap_elements(p5) | wrap_elements(p10)
