Load packages
## add 'developer' to packages to be installed from github
x <- c("data.table", "lubridate", "devtools", "maRce10/warbleR", "readxl", "ranger", "caret", "e1071", "pbapply", "viridis", "ggplot2", "kableExtra", "rlang", "Sim.DiffProc", "lme4", "soundgen", "markovchain", "igraph"#, "TraMineR", "spgs"
)
aa <- lapply(x, function(y) {
# get pakage name
pkg <- strsplit(y, "/")[[1]]
pkg <- pkg[length(pkg)]
# check if installed, if not then install
if (!pkg %in% installed.packages()[,"Package"]) {
if (grepl("/", y)) devtools::install_github(y, force = TRUE) else
install.packages(y)
}
# load package
a <- try(require(pkg, character.only = T), silent = T)
if (!a) remove.packages(pkg)
})
Functions and global parameters
warbleR_options(wl = 300, parallel = 1, bp = "frange", fast = TRUE, threshold = 15, ovlp = 20)
opts_knit$set(root.dir = "..")
# set evaluation false
opts_chunk$set( fig.width = 6, fig.height = 3, eval = TRUE, warning = FALSE, message = FALSE, tidy = TRUE)
clls <- readRDS("./data/processed/curated_extended_selection_table_inquiry_calls_2020_&_2021.RDS")
metadat <- read.csv("./data/processed/metadata_inquiry_calls_2020_&_2021.csv")
drs <- list.dirs("/home/m/Dropbox/Salidas de videos analizados/Salidas50FramesXBicho", full.names = FALSE)
drs <- drs[2:19]
videos <- sapply(strsplit(drs, "_"), "[[", 1)
videos <- as.numeric(gsub("T", "", videos))
years <- as.numeric(sapply(strsplit(drs, "_"), "[[", 3))
year.video <- paste(videos, years, sep = "-")
video_df <- data.frame(directories = drs, video = videos, year = years, year.video = year.video)
metadat$year.video <- paste(metadat$Video, metadat$year, sep = "-")
video_df$year.audio <- sapply(1:nrow(video_df), function(x){
metadat$year.audio[metadat$year.video == video_df$year.video[x]][1]
})
metadat$year.video <- paste(metadat$Video, metadat$year, sep = "-")
video_df$experiment <- sapply(1:nrow(video_df), function(x){
metadat$Experimento[metadat$year.video == video_df$year.video[x]][1]
})
video_df$flight.time <- sapply(1:nrow(video_df), function(x){
metadat$Tiempo.de.vuelo[metadat$year.video == video_df$year.video[x]][1]
})
# fix flight time
video_df$flight.time <- video_df$flight.time * 1440
video_df$n.calls <- sapply(1:nrow(video_df), function(x){
sum(clls$year.audio == video_df$year.audio[x])
})
att.cll <- attr(clls, "check.results")
# get call rate for the frist 2 minutes
video_df$n.calls.2.min <- sapply(1:nrow(video_df), function(x){
sum(clls$year.audio == video_df$year.audio[x] & att.cll$orig.start <= 120)
})
video_df$call.rate.2.min <- video_df$n.calls.2.min / video_df$flight.time
boxplot(call.rate.2.min ~ factor(flight.time), data = video_df)
fls <- list.files("/home/m/Dropbox/Salidas de videos analizados/Salidas50FramesXBicho", full.names = TRUE, recursive = TRUE)
fls <- grep("unpaired", fls, value = TRUE)
fls <- grep("volver", invert = TRUE, fls, value = TRUE)
video.3d.coords_l <- lapply(fls, function(x){
y <- read.csv(x)
y$video <- basename(x)
return(y)
})
names(video.3d.coords_l) <- basename(fls)
dists_by_video <- lapply(video.3d.coords_l, function(x){
coor_df <- data.frame(x = stack(x[, grep("x", names(x))])$values, y = stack(x[, grep("y", names(x))])$values, z = stack(x[, grep("z", names(x))])$values, ind = rep(1:((ncol(x) -1) / 3), each = nrow(x)), frame = 1:nrow(x))
dists <- sapply(unique(coor_df$frame), function(w){
W <- coor_df[coor_df$frame == w, ]
dst <- dist(W[, !names(coor_df) %in% c("ind", "frame")])
return(mean(dst)[1])
})
out <- data.frame(video = x$video[1], frame = unique(coor_df$frame), distance = dists, n.ind = length(unique(coor_df$ind)))
return(out)
})
dists_by_video_df <- do.call(rbind,
dists_by_video)
dists_by_video_df$year <- sapply(strsplit(dists_by_video_df$video, "_"), "[[", 3)
dists_by_video_df$year <- sapply(strsplit(dists_by_video_df$year, "-"), "[[", 1)
dists_by_video_df$video <- sapply(strsplit(dists_by_video_df$video, "_"), "[[", 1)
dists_by_video_df$video <- as.numeric(gsub("T", "", dists_by_video_df$video))
dists_by_video_df$year.video <- paste(dists_by_video_df$video, dists_by_video_df$year, sep = "-")
rownames(dists_by_video_df) <- 1:nrow(dists_by_video_df)
agg_dist <- aggregate(cbind(distance, n.ind) ~ year.video, data = dists_by_video_df, mean)
dist_callrate_df <- merge(agg_dist, video_df[, c("year.video", "call.rate.2.min")])
dist_callrate_df$n.ind.f <- factor(dist_callrate_df$n.ind, sort(unique(dist_callrate_df$n.ind)))
ggplot(dist_callrate_df, aes(x = call.rate.2.min, y = distance, size = n.ind, col = n.ind)) +
geom_point() +
scale_color_viridis_c(end = 0.8) +
labs( x = "Call rate (calls / s)", y = "Mean distance (m)", col = "Group\n size") +
xlim(c(0, 0.37)) +
geom_smooth(method = "lm") +
guides(size = "none") +
theme_classic(base_size = 30) +
theme(legend.position = c(0.94, 0.5),
legend.direction = "vertical")
mod1 <- lmer(distance ~ call.rate.2.min + (1|n.ind), dist_callrate_df, control = lmerControl(optimizer= "optimx", optCtrl = list(method="nlminb")))
summary(mod1)
## Linear mixed model fit by REML ['lmerMod']
## Formula: distance ~ call.rate.2.min + (1 | n.ind)
## Data: dist_callrate_df
## Control: lmerControl(optimizer = "optimx", optCtrl = list(method = "nlminb"))
##
## REML criterion at convergence: 28
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.44664 -0.72067 -0.00758 0.65053 1.74815
##
## Random effects:
## Groups Name Variance Std.Dev.
## n.ind (Intercept) 0.0000 0.0000
## Residual 0.2854 0.5342
## Number of obs: 18, groups: n.ind, 5
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 1.6586 0.1695 9.788
## call.rate.2.min 0.1132 0.6146 0.184
##
## Correlation of Fixed Effects:
## (Intr)
## cll.rt.2.mn -0.669
## optimizer (optimx) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular
mod2 <- lm(distance ~ call.rate.2.min + n.ind, data = dist_callrate_df)
summary(mod2)
##
## Call:
## lm(formula = distance ~ call.rate.2.min + n.ind, data = dist_callrate_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.67717 -0.37743 -0.03094 0.21990 1.08502
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.343371 0.392364 3.424 0.00377 **
## call.rate.2.min -0.006759 0.632964 -0.011 0.99162
## n.ind 0.082069 0.091986 0.892 0.38638
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5377 on 15 degrees of freedom
## Multiple R-squared: 0.0524, Adjusted R-squared: -0.07395
## F-statistic: 0.4147 on 2 and 15 DF, p-value: 0.6679
mod3 <- lm(distance ~ factor(n.ind), data = dist_callrate_df)
summary(mod3)
##
## Call:
## lm(formula = distance ~ factor(n.ind), data = dist_callrate_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.64747 -0.39570 -0.04792 0.22135 1.03888
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.55365 0.33217 4.677 0.000433 ***
## factor(n.ind)3 -0.03674 0.43942 -0.084 0.934640
## factor(n.ind)4 0.15820 0.46976 0.337 0.741673
## factor(n.ind)5 0.21060 0.43942 0.479 0.639712
## factor(n.ind)6 0.27389 0.43942 0.623 0.543862
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5753 on 13 degrees of freedom
## Multiple R-squared: 0.05961, Adjusted R-squared: -0.2297
## F-statistic: 0.206 on 4 and 13 DF, p-value: 0.9305