# package names
packages<-c("dplyr","readr","ggplot2", "here", "mcp", "knitr", "ezknitr", "loo")
# install any packages not previously installed
installed_packages<-packages %in% rownames(installed.packages())
if(any(installed_packages == FALSE)){
install.packages(packages[!installed_packages])
}
# load packages
invisible(lapply(packages, library, character.only = TRUE))
Import the first 10 days of net-squared displacement (NSD) data for trumpeter swan 3E.
df<-read_csv(here("DRUM_materials/data_files_for_DRUM/post_capture_swan.csv"))
Calculate average hourly NSD values
df<-df %>%
group_by(hours_since_capture) %>%
mutate(hr_nsd=mean(nsd)) %>%
distinct(hours_since_capture, hr_nsd)
Now the dataset is thinned to a size that will be easier to process. Try a model with just one intercept, representing no capture response, and two intercepts, indicating a capture response.
m0=list(hr_nsd~1)
m1=list(
hr_nsd~1,
~1)
f1<-mcp(m0, df[,c("hr_nsd","hours_since_capture")],
par_x="hours_since_capture")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 258
## Unobserved stochastic nodes: 2
## Total graph size: 1309
##
## Initializing model
f2<-mcp(m1, df[,c("hr_nsd","hours_since_capture")],
par_x="hours_since_capture")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 258
## Unobserved stochastic nodes: 4
## Total graph size: 3116
##
## Initializing model
Check models with cross-validation
f1$loo<-loo(f1)
f2$loo<-loo(f2)
loo_compare(f1$loo, f2$loo)
## elpd_diff se_diff
## model2 0.0 0.0
## model1 -215.1 62.6
The model with two intercepts out-performs the model with just one.
Plot the data.
plot(f2, q_fit=T)+
ggtitle("Trumpeter Swan flee and return post-capture")+
xlab("Hours since capture")+ylab("Net-squared displacement\n")+
theme_bw(base_size = 16)+
theme(plot.title=element_text(hjust=0.5))
sessionInfo()
## R version 4.0.2 (2020-06-22)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 19043)
##
## Matrix products: default
##
## locale:
## [1] LC_COLLATE=English_United States.1252
## [2] LC_CTYPE=English_United States.1252
## [3] LC_MONETARY=English_United States.1252
## [4] LC_NUMERIC=C
## [5] LC_TIME=English_United States.1252
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] ggplot2_3.3.2 loo_2.4.1 ezknitr_0.6 knitr_1.30 mcp_0.3.0.9000
## [6] here_0.1 readr_1.3.1 dplyr_1.0.2
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.1.0 xfun_0.26 bslib_0.3.1
## [4] purrr_0.3.4 lattice_0.20-41 colorspace_1.4-1
## [7] vctrs_0.3.4 generics_0.0.2 htmltools_0.5.2
## [10] yaml_2.2.1 rlang_0.4.10 jquerylib_0.1.4
## [13] pillar_1.4.6 glue_1.4.2 withr_2.3.0
## [16] HDInterval_0.2.2 distributional_0.2.2 plyr_1.8.6
## [19] matrixStats_0.57.0 lifecycle_0.2.0 stringr_1.4.0
## [22] munsell_0.5.0 gtable_0.3.0 codetools_0.2-16
## [25] coda_0.19-3 evaluate_0.14 labeling_0.3
## [28] forcats_0.5.0 fastmap_1.1.0 parallel_4.0.2
## [31] Rcpp_1.0.8 arrayhelpers_1.1-0 backports_1.1.10
## [34] scales_1.1.1 checkmate_2.0.0 jsonlite_1.7.1
## [37] farver_2.0.3 tidybayes_2.3.1 hms_0.5.3
## [40] digest_0.6.28 svUnit_1.0.3 stringi_1.5.3
## [43] grid_4.0.2 rprojroot_1.3-2 ggdist_2.4.0
## [46] tools_4.0.2 magrittr_2.0.1 sass_0.4.0
## [49] patchwork_1.1.1 tibble_3.0.3 tidyr_1.1.2
## [52] crayon_1.3.4 pkgconfig_2.0.3 ellipsis_0.3.1
## [55] rmarkdown_2.11 rstudioapi_0.11 R6_2.4.1
## [58] compiler_4.0.2