Echolocation plasticity

Author

Marcelo Araya-Salas

Published

June 19, 2023

1 Purpose

  • The first goal of this report

  • The second goal of this report

 

2 Report overview

 

3 Load packages

Code 
# knitr is require for creating html/pdf/word reports formatR is
# used for soft-wrapping code

# install/ load packages
sketchy::load_packages(packages = c("knitr", "formatR", "warbleR",
    "Rraven", "ohun", "readxl", "brms", "brmsish", "viridis", "ggplot2"))
Warning: replacing previous import 'brms::rstudent_t' by 'ggdist::rstudent_t'
when loading 'brmsish'
Warning: replacing previous import 'brms::pstudent_t' by 'ggdist::pstudent_t'
when loading 'brmsish'
Warning: replacing previous import 'brms::dstudent_t' by 'ggdist::dstudent_t'
when loading 'brmsish'
Warning: replacing previous import 'brms::qstudent_t' by 'ggdist::qstudent_t'
when loading 'brmsish'
Warning: replacing previous import 'brms::rhat' by 'posterior::rhat' when
loading 'brmsish'

4 Clipping recordings

Code 
anns <- imp_raven(path = "./data/raw/clip_annotations", warbler.format = TRUE,
    all.data = TRUE, name.from.file = TRUE, ext.case = "lower")


feature_reference(anns, units = c("s", "kHz"))

cut_sels(anns, path = "./data/raw/recordings", dest.path = "./data/raw/clips")


warbleR_options(wav.path = "./data/raw/clips")


full_spectrograms(dest.path = "./data/processed/spectrograms", flim = c(20,
    120), sxrow = 0.25, rows = 8, fast.spec = FALSE, ovlp = 90, X = NULL,
    horizontal = TRUE, parallel = 4, collevels = seq(-100, 0, 5),
    suffix = "0.25s")

5 Importing annotations

Code 
warbleR_options(wav.path = "./data/raw/clips")

anns <- imp_raven(path = "./data/processed/annotations", warbler.format = TRUE,
    all.data = TRUE, name.from.file = TRUE, ext.case = "lower")

anns$call_type <- ifelse(grepl("Table", anns$sound.files), "B", "A")


anns$sound.files <- gsub(".Table.1.selections 2", "", anns$sound.files)


cs <- check_sels(anns, fix.selec = TRUE)

sp_feat <- spectro_analysis(anns, ovlp = 70)
sp_feat$modindx <- NULL

sp_feat$pc1 <- prcomp(sp_feat[, -c(1:2)], scale. = TRUE)$x[, 1]


sp_feat$treatment <- ifelse(grepl("close", sp_feat$sound.files, fixed = TRUE),
    "closed", "open")

sp_feat$individual <- sapply(strsplit(sp_feat$sound.files, "_"), "[",
    1)

metadata <- read_excel("./data/raw/Bats_data.xlsx")

sp_feat$genus <- sapply(seq_len(nrow(sp_feat)), function(x) metadata$species[metadata$individual ==
    sp_feat$individual[x]][1])

write.csv(sp_feat, "./data/processed/acoustic_features.csv", row.names = FALSE)

6 Statistical analysis

Code 
sp_feat <- read.csv("./data/processed/acoustic_features.csv")


iter <- 10000
chains <- 4
priors <- c(prior(normal(0, 6), class = "b"), prior(cauchy(0, 4),
    class = "sd"))

pc1_mod <- brm(formula = pc1 ~ treatment + (1 | individual) + (1 |
    genus), data = sp_feat, prior = priors, chains = chains, iter = iter)


saveRDS(pc1_mod, "./data/processed/pc1_model.RDS")


freq_mod <- brm(formula = meanpeakf ~ treatment + (1 | individual) +
    (1 | genus), data = sp_feat, prior = priors, chains = chains,
    iter = iter)


saveRDS(freq_mod, "./data/processed/freq_model.RDS")


dur_mod <- brm(formula = duration ~ treatment + (1 | individual) +
    (1 | genus), data = sp_feat, prior = priors, chains = chains,
    iter = iter)


saveRDS(dur_mod, "./data/processed/duration_model.RDS")


bw_mod <- brm(formula = freq.IQR ~ treatment + (1 | individual) +
    (1 | genus), data = sp_feat, prior = priors, chains = chains,
    iter = iter)


saveRDS(bw_mod, "./data/processed/bandwidth_model.RDS")
Code 
source("~/Dropbox/R_package_testing/brmsish/R/extended_summary.R")
source("~/Dropbox/R_package_testing/brmsish/R/helpers.R")

my.viridis <- function(...) viridis(alpha = 0.5, begin = 0.3, end = 0.7,
    ...)

extended_summary(read.file = "./data/processed/pc1_model.RDS", n.posterior = 1000,
    fill = "orange3", trace.palette = my.viridis, remove.intercepts = FALSE,
    highlight = TRUE)

6.1 pc1_model

priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 0, 3.2) sd-cauchy(0, 4) sigma-student_t(3, 0, 3.2) pc1 ~ treatment + (1 | individual) + (1 | genus) 10000 4 1 5000 302 0 2039.008 1117.255 1653856741
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
b_Intercept 1.309 -1.278 3.856 1.002 2039.008 1117.255
b_treatmentopen -1.558 -2.857 -0.227 1 11808.204 6009.440

Code 
extended_summary(read.file = "./data/processed/freq_model.RDS", n.posterior = 1000,
    fill = "orange3", trace.palette = my.viridis, remove.intercepts = FALSE,
    highlight = TRUE)

6.2 freq_model

priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 77.7, 12) sd-cauchy(0, 4) sigma-student_t(3, 0, 12) meanpeakf ~ treatment + (1 | individual) + (1 | genus) 5000 4 1 2500 40 0 2885.746 2385.547 2070648850
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
b_Intercept 81.929 75.105 88.952 1.001 2885.746 2385.547
b_treatmentopen -5.121 -9.528 -0.610 1 8028.631 5509.010

Code 
extended_summary(read.file = "./data/processed/duration_model.RDS",
    n.posterior = 1000, fill = "orange3", trace.palette = my.viridis,
    remove.intercepts = FALSE, highlight = TRUE)

6.3 duration_model

priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 0, 2.5) sd-cauchy(0, 4) sigma-student_t(3, 0, 2.5) duration ~ treatment + (1 | individual) + (1 | genus) 5000 4 1 2500 35 0 5721.414 5681.433 1682270005
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
b_Intercept 0.001 0.001 0.002 1.001 5721.414 5681.433
b_treatmentopen 0.000 0.000 0.000 1 10698.867 8276.534

Code 
extended_summary(read.file = "./data/processed/bandwidth_model.RDS",
    n.posterior = 1000, fill = "orange3", trace.palette = my.viridis,
    remove.intercepts = FALSE, highlight = TRUE)

6.4 bandwidth_model

priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 9.3, 2.9) sd-cauchy(0, 4) sigma-student_t(3, 0, 2.9) freq.IQR ~ treatment + (1 | individual) + (1 | genus) 5000 4 1 2500 133 0 402.545 119.691 1101101355
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
b_Intercept 11.841 8.386 16.317 1.016 402.545 119.691
b_treatmentopen -1.645 -3.592 0.294 1.002 2909.844 4749.666

Code 
sp_feat <- read.csv("./data/processed/acoustic_features.csv")

ggplot(sp_feat, aes(x = treatment, y = pc1)) + geom_boxplot(fill = viridis(10,
    alpha = 0.3)[4]) + labs(x = "Treatment", y = "PC1")

Code 
ggplot(sp_feat, aes(x = treatment, y = meanpeakf)) + geom_boxplot(fill = viridis(10,
    alpha = 0.3)[4]) + labs(x = "Treatment", y = "Peak frequency (kHz)")

Code 
ggplot(sp_feat, aes(x = treatment, y = duration)) + geom_boxplot(fill = viridis(10,
    alpha = 0.3)[4]) + labs(x = "Treatment", y = "Duration (s)")

Code 
ggplot(sp_feat, aes(x = treatment, y = freq.IQR)) + geom_boxplot(fill = viridis(10,
    alpha = 0.3)[4]) + labs(x = "Treatment", y = "Bandwidth (kHz)")

Takeaways

Only peak frequency differ between open and closed treatments: individuals increase their frequency when flying in a closed (cluttered) environment

 

 

Session information

R version 4.2.2 Patched (2022-11-10 r83330)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 20.04.5 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0

locale:
 [1] LC_CTYPE=es_ES.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=es_CR.UTF-8        LC_COLLATE=es_ES.UTF-8    
 [5] LC_MONETARY=es_CR.UTF-8    LC_MESSAGES=es_ES.UTF-8   
 [7] LC_PAPER=es_CR.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=es_CR.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] ggplot2_3.4.2      viridis_0.6.3      viridisLite_0.4.2  brmsish_1.0.0     
 [5] brms_2.18.0        Rcpp_1.0.10        readxl_1.4.1       ohun_0.1.1        
 [9] Rraven_1.0.14      warbleR_1.1.28     NatureSounds_1.0.4 seewave_2.2.0     
[13] tuneR_1.4.4        formatR_1.12       knitr_1.42        

loaded via a namespace (and not attached):
  [1] backports_1.4.1      systemfonts_1.0.4    plyr_1.8.7          
  [4] igraph_1.3.5         crosstalk_1.2.0      rstantools_2.2.0    
  [7] inline_0.3.19        digest_0.6.31        htmltools_0.5.5     
 [10] fansi_1.0.4          magrittr_2.0.3       checkmate_2.1.0     
 [13] remotes_2.4.2        RcppParallel_5.1.5   matrixStats_0.62.0  
 [16] xts_0.12.2           svglite_2.1.0        prettyunits_1.1.1   
 [19] colorspace_2.1-0     signal_0.7-7         rvest_1.0.3         
 [22] ggdist_3.2.1         xfun_0.39            dplyr_1.1.0         
 [25] callr_3.7.3          crayon_1.5.2         RCurl_1.98-1.12     
 [28] jsonlite_1.8.4       zoo_1.8-11           ape_5.7-1           
 [31] glue_1.6.2           kableExtra_1.3.4     gtable_0.3.3        
 [34] emmeans_1.8.6        webshot_0.5.4        distributional_0.3.1
 [37] pkgbuild_1.3.1       rstan_2.21.7         abind_1.4-5         
 [40] scales_1.2.1         mvtnorm_1.1-3        DBI_1.1.3           
 [43] xaringanExtra_0.7.0  miniUI_0.1.1.1       dtw_1.23-1          
 [46] xtable_1.8-4         units_0.8-1          proxy_0.4-27        
 [49] stats4_4.2.2         StanHeaders_2.21.0-7 DT_0.26             
 [52] htmlwidgets_1.5.4    httr_1.4.6           threejs_0.3.3       
 [55] posterior_1.3.1      ellipsis_0.3.2       pkgconfig_2.0.3     
 [58] loo_2.5.1            farver_2.1.1         utf8_1.2.3          
 [61] labeling_0.4.2       tidyselect_1.2.0     rlang_1.1.1         
 [64] reshape2_1.4.4       later_1.3.1          munsell_0.5.0       
 [67] cellranger_1.1.0     tools_4.2.2          cli_3.6.1           
 [70] generics_0.1.3       ggridges_0.5.4       evaluate_0.21       
 [73] stringr_1.5.0        fastmap_1.1.1        yaml_2.3.7          
 [76] processx_3.8.1       packrat_0.8.1        pbapply_1.7-0       
 [79] nlme_3.1-162         mime_0.12            xml2_1.3.4          
 [82] brio_1.1.3           compiler_4.2.2       bayesplot_1.9.0     
 [85] shinythemes_1.2.0    rstudioapi_0.14      e1071_1.7-12        
 [88] testthat_3.1.8       sketchy_1.0.2        tibble_3.2.1        
 [91] stringi_1.7.12       ps_1.7.5             Brobdingnag_1.2-9   
 [94] lattice_0.20-45      Matrix_1.5-1         classInt_0.4-8      
 [97] markdown_1.3         shinyjs_2.1.0        fftw_1.0-7          
[100] tensorA_0.36.2       vctrs_0.6.2          pillar_1.9.0        
[103] lifecycle_1.0.3      bridgesampling_1.1-2 estimability_1.4.1  
[106] cowplot_1.1.1        bitops_1.0-7         httpuv_1.6.6        
[109] R6_2.5.1             promises_1.2.0.1     KernSmooth_2.23-20  
[112] gridExtra_2.3        codetools_0.2-19     colourpicker_1.2.0  
[115] MASS_7.3-58.2        gtools_3.9.3         rjson_0.2.21        
[118] withr_2.5.0          shinystan_2.6.0      parallel_4.2.2      
[121] grid_4.2.2           coda_0.19-4          class_7.3-21        
[124] rmarkdown_2.21       sf_1.0-9             shiny_1.7.3         
[127] base64enc_0.1-3      dygraphs_1.1.1.6