library(knitr) #for R Markdown
library(MASS) #for negative binomial glm (glm.nb)
library(lme4) #for mixed models
library(emmeans) #for posthoc
library(car) #for Anova
library(survival) #for survival analysis
library(coxme)
library(rptR) #for repeatability analysis
library(MuMIn) #for model selection (dredge)
library(ggfortify) #for plotting survival analysis
library(ggsignif) #for labeling significance in ggplots
library(GGally) #for correlation matrix
library(tidyverse) #for data processing, put last to avoid function masking
data_raw <- read.csv("data_maze_040326.csv",
# to avoid reading errors
fileEncoding="UTF-8-BOM", na.strings = "")
data <- data_raw %>%
#drop trailing NA reads
filter(!if_all(everything(), is.na)) %>%
#correct variable types
mutate_at(c("ID","sex","assayer","trial","shelter_YN","shelter_quadrant"), as.factor) %>%
mutate_at(c("found_date","morph_date","assay_date"), lubridate::mdy) %>%
#derive variables
mutate(
#life history variables
lat_meta = as.numeric(morph_date - found_date),
age = as.numeric(assay_date - morph_date),
#performance variables
lat_shelter_from_movement = as.numeric(lat_shelter - lat_move),
lines_rate = lines_crossed/lat_trial,
grids_rate = grids_explored/lat_trial,
#for survival analysis and plotting
shelter_10 = recode(shelter_YN, "Y" = 1, "N" = 0),
move_10 = ifelse(lat_move==lat_trial, 0, 1)
)
Make an aggregrate version of the dataset by frog ID
data_ind <- data %>%
group_by(ID) %>%
summarise(lat_move_avg = mean(lat_move),
lat_shelter_avg = mean(lat_shelter),
lat_shelter_mov_avg = mean(lat_shelter_from_movement),
lat_shelter_improv = lat_shelter[trial == 1] - lat_shelter[trial == 5],
lines_rate_avg = mean(lines_rate),
grids_rate_avg = mean(grids_rate),
prop_complete = sum(shelter_YN == "Y")/n(),
found_date = mean(found_date),
morph_date = mean(morph_date),
age = mean(age),
lat_meta = mean(lat_meta)
)
rpt{rptR}: Repeatability EstimationP_permut: p-value from permutation test (more
robust?)LRT_P: p-value from likelihood ratio testAll three personality axes are repeatable!
rpt_result <- rpt(lat_move ~ (1 | ID),
grname = "ID",
data = data,
datatype = "Gaussian",
nboot = 1000,
npermut = 1000)
## Bootstrap Progress:
## Permutation Progress for ID :
summary(rpt_result)$rpt %>% kable(digits = 3)
|
rpt_result <- rpt(lines_rate ~ (1 | ID),
grname = "ID",
data = data,
datatype = "Gaussian",
nboot = 1000,
npermut = 1000)
## Bootstrap Progress:
## Permutation Progress for ID :
summary(rpt_result)$rpt %>% kable(digits = 3)
|
rpt_result <- rpt(grids_rate ~ (1 | ID),
grname = "ID",
data = data,
datatype = "Gaussian",
nboot = 1000,
npermut = 1000)
## Bootstrap Progress:
## Permutation Progress for ID :
summary(rpt_result)$rpt %>% kable(digits = 3)
|
mod_boldness_activity <- lmer(lines_rate ~ lat_move + (1|ID),
data = data)
summary(mod_boldness_activity)$coefficients
## Estimate Std. Error t value
## (Intercept) 2.224513e-02 2.225872e-03 9.9938964
## lat_move -1.150276e-05 1.265838e-05 -0.9087068
Anova(mod_boldness_activity)
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: lines_rate
## Chisq Df Pr(>Chisq)
## lat_move 0.8257 1 0.3635
mod_boldness_exploration <- lmer(grids_rate ~ lat_move + (1|ID),
data = data)
summary(mod_boldness_exploration)$coefficients
## Estimate Std. Error t value
## (Intercept) 7.529567e-03 6.981217e-04 10.785465
## lat_move 5.257119e-06 3.647541e-06 1.441277
Anova(mod_boldness_exploration)
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: grids_rate
## Chisq Df Pr(>Chisq)
## lat_move 2.0773 1 0.1495
mod_activity_exploration <- lmer(grids_rate ~ lines_rate + (1|ID),
data = data)
summary(mod_activity_exploration)$coefficients
## Estimate Std. Error t value
## (Intercept) 0.005546651 0.0008851838 6.266101
## lines_rate 0.110757968 0.0354561117 3.123805
Anova(mod_activity_exploration)
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: grids_rate
## Chisq Df Pr(>Chisq)
## lines_rate 9.7582 1 0.001785 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
fig_boldness_activity <-
ggplot(data, aes(x = lat_move, y = lines_rate)) +
geom_point(color = "black", size = 2) +
geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("Latency to move (s)") +
ylab("# Lines crossed / Trial time (s)") +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
fig_boldness_exploration <-
ggplot(data, aes(x = lat_move, y = grids_rate)) +
geom_point(color = "black", size = 2) +
geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("Latency to move (s)") +
ylab("# Grids explored / Trial time (s)") +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
fig_activity_exploration <-
ggplot(data, aes(x = lines_rate, y = grids_rate)) +
geom_point(color = "black", size = 2) +
geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("# Lines crossed / Trial time (s)") +
ylab("# Grids explored / Trial time (s)") +
theme_bw(base_size = 15) +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
egg::ggarrange(fig_boldness_activity, fig_boldness_exploration, fig_activity_exploration,
nrow = 1,
labels = c("A","B","C"))
fig_boldness_activity_ind <-
ggplot(data_ind, aes(x = lat_move_avg, y = lines_rate_avg)) +
# for overlaying trial level data - doesn't look good on the plot though
# geom_point(data = data, aes(x = lat_move, y = lines_rate),
# color = "grey80", size = 1) +
geom_point(color = "black", size = 2) +
geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("Latency to move (s)") +
ylab("# Lines crossed / Trial time (s)") +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
fig_boldness_exploration_ind <-
ggplot(data_ind, aes(x = lat_move_avg, y = grids_rate_avg)) +
geom_point(color = "black", size = 2) +
geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("Latency to move (s)") +
ylab("# Grids explored / Trial time (s)") +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
fig_activity_exploration_ind <-
ggplot(data_ind, aes(x = lines_rate_avg, y = grids_rate_avg)) +
geom_point(color = "black", size = 2) +
geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("# Lines crossed / Trial time (s)") +
ylab("# Grids explored / Trial time (s)") +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
egg::ggarrange(fig_boldness_activity_ind, fig_boldness_exploration_ind, fig_activity_exploration_ind,
nrow = 1,
labels = c("A","B","C"))
mod<- lm(prop_complete ~ lat_meta, data = data_ind)
summary(mod)$coefficients
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.51314286 1.10907807 -0.4626751 0.6535013
## lat_meta 0.01314286 0.01866657 0.7040853 0.4974475
mod<- lm(prop_complete ~ age, data = data_ind)
summary(mod)$coefficients
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.1186611406 0.613392068 0.1934507 0.8504788
## age 0.0003952963 0.001630547 0.2424317 0.8133454
fig_meta_comp <-
ggplot(data_ind ,
aes(x = lat_meta, y = prop_complete)) +
geom_point(color = "black", size = 2) +
#geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("Latency to Metamorphose (days)") +
ylab("Percentage of Success (%)") +
scale_y_continuous(labels = scales::percent) +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
fig_age_comp <-
ggplot(data_ind,
aes(x = age, y = prop_complete)) +
geom_point(color = "black", size = 2) +
# geom_smooth(method = lm, alpha = 0.2, color = "SteelBlue", fill = "SteelBlue")+
xlab("Age (days)") +
ylab("Percentage of Success (%)") +
scale_y_continuous(labels = scales::percent) +
theme_bw(base_size = 15)+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
aspect.ratio=1)
egg::ggarrange(fig_meta_comp, fig_age_comp,
nrow = 1,
labels = c("A","B"))
No clear sign of learning! They overall did best at finding shelter in the first trial.
data_plot <-
survfit(Surv(lat_shelter, shelter_10) ~ trial,
data = data) %>% #filter out the accidental 6th trials
fortify(surv.connect = TRUE) #make into a dataframe
fig_shelter_lat <-
ggplot(data_plot, aes(time, 1-surv, color = strata)) +
#survival plot
geom_step(linewidth = 1) +
# axes and legends
labs(x = "Time since start of trial (s)", y = "Proportion Frogs in Shelter") +
scale_color_manual(values = rev(hcl.colors(7, palette = "BuPu")[1:5]),
name = "Trial number") +
scale_y_continuous(expand = c(0,0), limit = c(0,0.61),
labels = scales::percent) +
# adjust element themes
theme_classic(base_size = 15)
fig_shelter_lat