GLMMs
Global model
coor_res <- read.csv("HOWR coord by focal bird.csv")
coor_res <- coor_res[complete.cases(coor_res),]
# remove the ones with less than 20 obs for each focal and no focal
# coor_res <- coor_res[coor_res$n.focal > 20 & coor_res$n.no.focal > 20, ]
# normalize
pre.dat <- preProcess(coor_res[, sapply(coor_res, is.numeric)], method = c("BoxCox", "scale", "center"))
coor_res[, sapply(coor_res, is.numeric)] <- predict(pre.dat, coor_res[, sapply(coor_res, is.numeric)])
nrow(coor_res)## [1] 259glob_mod <- lmer(coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr), data = coor_res, REML = FALSE)## Warning: 'rBind' is deprecated.
## Since R version 3.2.0, base's rbind() should work fine with S4 objectsoptions(na.action = "na.fail")
drg <- dredge(global.model = glob_mod, extra = alist(AIC))## Fixed term is "(Intercept)"drg## Global model call: lmer(formula = coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr),
## data = coor_res, REML = FALSE)
## ---
## Model selection table
## (Intrc) adjtf noise trt AIC df logLik AICc delta weight
## 1 0.004551 738.9 3 -366.452 739.0 0.00 0.502
## 3 0.002850 -0.04460 740.5 4 -366.265 740.7 1.69 0.216
## 2 0.021440 + 740.7 4 -366.370 740.9 1.90 0.194
## 4 0.023510 + -0.04880 742.3 5 -366.144 742.5 3.53 0.086
## 5 0.256600 + 749.6 11 -363.798 750.7 11.67 0.001
## 7 0.261000 -0.04853 + 751.1 12 -363.573 752.4 13.42 0.001
## 6 0.278600 + + 751.4 12 -363.695 752.7 13.66 0.001
## 8 0.288000 + -0.05272 + 752.9 13 -363.425 754.3 15.34 0.000
## Models ranked by AICc(x)
## Random terms (all models):
## '1 | ntwrk.yr'glob_mod <- lmerTest::lmer(coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr), data = coor_res, REML = FALSE)
# summary(glob_mod)
lsm <- lmerTest::ls_means(glob_mod)
lsm## Least Squares Means table:
##
## Estimate Std. Error df t value lower
## trtduring noise 0.2506683 0.1825427 215.8 1.3732 -0.1091265
## trtduring song -0.0077655 0.1918300 215.2 -0.0405 -0.3858721
## trtduring song noise -0.1768622 0.1804027 217.2 -0.9804 -0.5324265
## trtpost noise -0.1136169 0.1890730 210.3 -0.6009 -0.4863385
## trtpost song 0.0338365 0.1773105 192.7 0.1908 -0.3158824
## trtpost song noise 0.0705669 0.1900262 219.9 0.3714 -0.3039388
## trtpre noise -0.0495372 0.2075527 224.8 -0.2387 -0.4585353
## trtpre song -0.2356901 0.1901922 220.0 -1.2392 -0.6105222
## trtpre song noise 0.0977218 0.1830238 211.6 0.5339 -0.2630620
## adjtofocalno 0.0228425 0.0809950 22.0 0.2820 -0.1451185
## adjtofocalyes -0.0518821 0.1150480 42.1 -0.4510 -0.2840346
## upper Pr(>|t|)
## trtduring noise 0.6104630 0.1711
## trtduring song 0.3703412 0.9677
## trtduring song noise 0.1787022 0.3280
## trtpost noise 0.2591046 0.5485
## trtpost song 0.3835555 0.8489
## trtpost song noise 0.4450726 0.7107
## trtpre noise 0.3594609 0.8116
## trtpre song 0.1391421 0.2166
## trtpre song noise 0.4585055 0.5940
## adjtofocalno 0.1908035 0.7806
## adjtofocalyes 0.1802703 0.6543
##
## Confidence level: 95%
## Degrees of freedom method: Satterthwaite# lmerTest::show_tests(lsm)Test model assumptions
qqnorm(resid(glob_mod), main="normal qq-plot, residuals", col = topo.colors(1), pch = 20, cex = 2)
qqline(resid(glob_mod), col = adjustcolor("red", alpha.f = 0.7))
sjstats::overdisp(glob_mod)##
## # Overdispersion test
##
## dispersion ratio = 0.9898
## Pearson's Chi-Squared = 243.4838
## p-value = 0.5334## No overdispersion detected.song models
coor_res <- read.csv("HOWR coord by focal bird.csv")
coor_res <- coor_res[complete.cases(coor_res),]
# remove the ones with less than 20 obs for each focal and no focal
# coor_res <- coor_res[coor_res$n.focal > 20 & coor_res$n.no.focal > 20, ]
coor_res <- coor_res[grep("song$", coor_res$trt), ]
# normalize
pre.dat <- preProcess(coor_res[, sapply(coor_res, is.numeric)], method = c("BoxCox", "scale", "center"))
coor_res[, sapply(coor_res, is.numeric)] <- predict(pre.dat, coor_res[, sapply(coor_res, is.numeric)])
nrow(coor_res)## [1] 89glob_mod <- lmer(coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr), data = coor_res, REML = FALSE)
options(na.action = "na.fail")
drg <- dredge(global.model = glob_mod, extra = alist(AIC))## Fixed term is "(Intercept)"drg## Global model call: lmer(formula = coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr),
## data = coor_res, REML = FALSE)
## ---
## Model selection table
## (Intrc) adjtf noise trt AIC df logLik AICc delta weight
## 1 8.732e-18 257.6 3 -125.783 257.8 0.00 0.290
## 5 7.892e-02 + 257.8 5 -123.905 258.5 0.68 0.206
## 2 6.086e-02 + 258.7 4 -125.329 259.1 1.29 0.153
## 6 1.270e-01 + + 258.9 6 -123.441 259.9 2.06 0.104
## 3 2.115e-14 -0.01557 259.5 4 -125.772 260.0 2.17 0.098
## 7 7.594e-02 0.02049 + 259.8 6 -123.886 260.8 2.95 0.066
## 4 6.108e-02 + -0.01765 260.6 5 -125.315 261.4 3.50 0.050
## 8 1.240e-01 + 0.01954 + 260.8 7 -123.424 262.2 4.38 0.032
## Models ranked by AICc(x)
## Random terms (all models):
## '1 | ntwrk.yr'glob_mod <- lmerTest::lmer(coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr), data = coor_res, REML = FALSE)
# summary(glob_mod)
lsm <- lmerTest::ls_means(glob_mod)
lsm## Least Squares Means table:
##
## Estimate Std. Error df t value lower upper
## trtduring song 0.012123 0.195214 89 0.0621 -0.375764 0.400009
## trtpost song 0.142316 0.177500 89 0.8018 -0.210373 0.495005
## trtpre song -0.343745 0.190792 89 -1.8017 -0.722844 0.035355
## adjtofocalno 0.048825 0.120280 89 0.4059 -0.190169 0.287818
## adjtofocalyes -0.175029 0.198758 89 -0.8806 -0.569956 0.219898
## Pr(>|t|)
## trtduring song 0.95062
## trtpost song 0.42482
## trtpre song 0.07498 .
## adjtofocalno 0.68577
## adjtofocalyes 0.38090
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Confidence level: 95%
## Degrees of freedom method: Satterthwaite# lmerTest::show_tests(lsm)Test model assumptions
qqnorm(resid(glob_mod), main="normal qq-plot, residuals", col = topo.colors(1), pch = 20, cex = 2)
qqline(resid(glob_mod), col = adjustcolor("red", alpha.f = 0.7))
sjstats::overdisp(glob_mod)##
## # Overdispersion test
##
## dispersion ratio = 1.0178
## Pearson's Chi-Squared = 83.4570
## p-value = 0.4344## No overdispersion detected.noise models
coor_res <- read.csv("HOWR coord by focal bird.csv")
coor_res <- coor_res[complete.cases(coor_res),]
# remove the ones with less than 20 obs for each focal and no focal
# coor_res <- coor_res[coor_res$n.focal > 20 & coor_res$n.no.focal > 20, ]
coor_res <- coor_res[grep("noise$", coor_res$trt), ]
# normalize
pre.dat <- preProcess(coor_res[, sapply(coor_res, is.numeric)], method = c("BoxCox", "scale", "center"))
coor_res[, sapply(coor_res, is.numeric)] <- predict(pre.dat, coor_res[, sapply(coor_res, is.numeric)])
nrow(coor_res)## [1] 170glob_mod <- lmer(coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr), data = coor_res, REML = FALSE)
options(na.action = "na.fail")
drg <- dredge(global.model = glob_mod, extra = alist(AIC))## Fixed term is "(Intercept)"drg## Global model call: lmer(formula = coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr),
## data = coor_res, REML = FALSE)
## ---
## Model selection table
## (Intrc) adjtf noise trt AIC df logLik AICc delta weight
## 1 0.0029630 487.0 3 -240.525 487.2 0.00 0.472
## 3 0.0009443 -0.07587 488.2 4 -240.109 488.5 1.27 0.251
## 2 0.0128900 + 489.0 4 -240.509 489.3 2.07 0.168
## 4 0.0178100 + -0.07828 490.1 5 -240.064 490.5 3.30 0.091
## 5 0.2041000 + 494.2 8 -239.120 495.1 7.94 0.009
## 7 0.2110000 -0.08665 + 495.2 9 -238.577 496.3 9.09 0.005
## 6 0.2130000 + + 496.2 9 -239.109 497.3 10.15 0.003
## 8 0.2281000 + -0.08879 + 497.1 10 -238.537 498.5 11.26 0.002
## Models ranked by AICc(x)
## Random terms (all models):
## '1 | ntwrk.yr'glob_mod <- lmerTest::lmer(coor.score ~ trt + adjtofocal + noise + (1 | ntwrk.yr), data = coor_res, REML = FALSE)Test model assumptions
qqnorm(resid(glob_mod), main="normal qq-plot, residuals", col = topo.colors(1), pch = 20, cex = 2)
qqline(resid(glob_mod), col = adjustcolor("red", alpha.f = 0.7))
sjstats::overdisp(glob_mod)##
## # Overdispersion test
##
## dispersion ratio = 1.0004
## Pearson's Chi-Squared = 160.0657
## p-value = 0.4837## No overdispersion detected.Session information
## R version 3.4.4 (2018-03-15)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 16.04.5 LTS
##
## Matrix products: default
## BLAS: /usr/lib/openblas-base/libblas.so.3
## LAPACK: /usr/lib/libopenblasp-r0.2.18.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] caret_6.0-78 lattice_0.20-35 MuMIn_1.40.4
## [4] lme4_1.1-15 Matrix_1.2-14 readxl_1.0.0
## [7] pbapply_1.3-4 knitr_1.20 kableExtra_0.9.0
## [10] ggplot2_2.2.1 warbleR_1.1.15 NatureSounds_1.0.0
## [13] seewave_2.1.0 tuneR_1.3.3 maps_3.3.0
##
## loaded via a namespace (and not attached):
## [1] TH.data_1.0-8 minqa_1.2.4 colorspace_1.3-2
## [4] rjson_0.2.20 class_7.3-14 modeltools_0.2-21
## [7] ggridges_0.5.0 sjlabelled_1.0.10 rprojroot_1.3-2
## [10] snakecase_0.9.1 estimability_1.3 rstudioapi_0.7
## [13] proxy_0.4-22 glmmTMB_0.2.0 DRR_0.0.3
## [16] Deriv_3.8.4 prodlim_1.6.1 mvtnorm_1.0-8
## [19] lubridate_1.7.1 coin_1.2-2 xml2_1.2.0
## [22] codetools_0.2-15 splines_3.4.4 mnormt_1.5-5
## [25] robustbase_0.92-8 sjmisc_2.7.1 RcppRoll_0.2.2
## [28] bayesplot_1.5.0 nloptr_1.0.4 broom_0.4.3
## [31] ddalpha_1.3.1 kernlab_0.9-25 sfsmisc_1.1-1
## [34] readr_1.1.1 compiler_3.4.4 httr_1.3.1
## [37] emmeans_1.1.3 sjstats_0.14.3 backports_1.1.2
## [40] assertthat_0.2.0 lazyeval_0.2.1 cli_1.0.0
## [43] htmltools_0.3.6 tools_3.4.4 bindrcpp_0.2.2
## [46] lmerTest_3.0-1 coda_0.19-1 gtable_0.2.0
## [49] glue_1.3.0 reshape2_1.4.3 dplyr_0.7.6
## [52] Rcpp_0.12.18 cellranger_1.1.0 nlme_3.1-137
## [55] iterators_1.0.10 lmtest_0.9-36 psych_1.7.8
## [58] timeDate_3042.101 gower_0.1.2 stringr_1.3.1
## [61] rvest_0.3.2 stringdist_0.9.4.6 DEoptimR_1.0-8
## [64] zoo_1.8-3 MASS_7.3-50 scales_0.5.0
## [67] ipred_0.9-6 hms_0.4.0 sandwich_2.4-0
## [70] pwr_1.2-2 TMB_1.7.13 Sim.DiffProc_4.0
## [73] yaml_2.1.19 rpart_4.1-13 stringi_1.2.2
## [76] foreach_1.4.4 e1071_1.6-8 lava_1.6
## [79] soundgen_1.2.1 rlang_0.2.1 pkgconfig_2.0.1
## [82] dtw_1.20-1 bitops_1.0-6 rgl_0.95.1441
## [85] pracma_2.1.4 evaluate_0.10.1 prediction_0.3.2
## [88] purrr_0.2.4 bindr_0.1.1 recipes_0.1.2
## [91] CVST_0.2-1 tidyselect_0.2.4 plyr_1.8.4
## [94] magrittr_1.5 R6_2.2.2 fftw_1.0-4
## [97] multcomp_1.4-8 dimRed_0.1.0 haven_1.1.1
## [100] pillar_1.2.3 foreign_0.8-70 withr_2.1.1
## [103] survival_2.42-6 scatterplot3d_0.3-40 RCurl_1.95-4.11
## [106] nnet_7.3-12 tibble_1.4.2 modelr_0.1.1
## [109] crayon_1.3.4 rmarkdown_1.9 jpeg_0.1-8
## [112] grid_3.4.4 data.table_1.10.4-3 forcats_0.3.0
## [115] ModelMetrics_1.1.0 digest_0.6.15 xtable_1.8-2
## [118] tidyr_0.7.2 numDeriv_2016.8-1 signal_0.7-6
## [121] stats4_3.4.4 munsell_0.4.3 viridisLite_0.3.0