Description
Asymptotic N-Mixture Model Fitting.
Description
For fitting N-mixture models using either FFT or asymptotic approaches. FFT N-mixture models extend the work of Cowen et al. (2017) <doi:10.1111/biom.12701>. Asymptotic N-mixture models extend the work of Dail and Madsen (2011) <doi:10.1111/j.1541-0420.2010.01465.x>, to consider asymptotic solutions to the open population N-mixture models. The FFT models are derived and described in "Parker, M.R.P., Elliott, L., Cowen, L.L.E. (2022). Computational efficiency and precision for replicated-count and batch-marked hidden population models [Manuscript in preparation]. Department of Statistics and Actuarial Sciences, Simon Fraser University.". The asymptotic models are derived and described in: "Parker, M.R.P., Elliott, L., Cowen, L.L.E., Cao, J. (2022). Fast asymptotic solutions for N-mixtures on large populations [Manuscript in preparation]. Department of Statistics and Actuarial Sciences, Simon Fraser University.".
README.md
quickNmix
The goal of quickNmix is to aid in the fitting of asymptotic N-mixture models, which are computed significantly faster than their canonical counterpart when population sizes are large.
Installation
You can install the released version of quickNmix from CRAN with:
install.packages("quickNmix")
And the development version from GitHub with:
# install.packages("devtools")
devtools::install_github("mrparker909/quickNmix")
Example 1
This is a basic example which shows how to fit a model with site varying lambda, and time varying pdet:
library(quickNmix)
tictoc::tic()
nit = anmu[1:2,1:5] # ancient murrelet chick counts
mod = fitNmix(nit=nit,
K=400, # upper bound on population size
l_s_c=list(c(0,1)), # lambda site covariate
p_t_c=list(c(0,1,1,1,1)),
control=list(reltol=1e-5))
#> Warning: executing %dopar% sequentially: no parallel backend registered
tictoc::toc()
#> 564.72 sec elapsed
model AIC value:
mod$model_results$AIC
#> [1] 77.5027
lambda estimates for each site:
mod$model_results$estimate_matrices$lambda
#> [,1]
#> [1,] 32.3202
#> [2,] 135.3340
gamma estimates for each site and time:
mod$model_results$estimate_matrices$gamma
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 25.84412 25.84412 25.84412 25.84412 25.84412
#> [2,] 25.84412 25.84412 25.84412 25.84412 25.84412
omega estimates for each site and time:
mod$model_results$estimate_matrices$omega
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.9572106 0.9572106 0.9572106 0.9572106 0.9572106
#> [2,] 0.9572106 0.9572106 0.9572106 0.9572106 0.9572106
pdet estimates for each site and time:
mod$model_results$estimate_matrices$pdet
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.9901256 0.4421962 0.4421962 0.4421962 0.4421962
#> [2,] 0.9901256 0.4421962 0.4421962 0.4421962 0.4421962
Example 2
This is a basic example which shows how you can use multiple cores to compute:
library(quickNmix)
# uses library doParallel, here we use 2 cores
doParallel::registerDoParallel(cores = 2)
tictoc::tic()
nit = anmu[c(2,5),1:12] # ancient murrelet chick counts
mod2 = fitNmix(nit=nit, K=400)
tictoc::toc()
#> 2908.69 sec elapsed
mod2
#> $optim_results
#> $optim_results$par
#> B_l_0 B_g_0 B_o_0 B_p_0
#> 4.995500 4.279050 -2.047628 3.029386
#>
#> $optim_results$value
#> [1] 107.5519
#>
#> $optim_results$counts
#> function gradient
#> 231 100
#>
#> $optim_results$convergence
#> [1] 1
#>
#> $optim_results$message
#> NULL
#>
#>
#> $model_results
#> $model_results$NLL
#> [1] 107.5519
#>
#> $model_results$AIC
#> [1] 223.1037
#>
#> $model_results$estimate_matrices
#> $model_results$estimate_matrices$lambda
#> [,1]
#> [1,] 147.7468
#> [2,] 147.7468
#>
#> $model_results$estimate_matrices$gamma
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 72.17185 72.17185 72.17185 72.17185 72.17185 72.17185 72.17185 72.17185
#> [2,] 72.17185 72.17185 72.17185 72.17185 72.17185 72.17185 72.17185 72.17185
#> [,9] [,10] [,11] [,12]
#> [1,] 72.17185 72.17185 72.17185 72.17185
#> [2,] 72.17185 72.17185 72.17185 72.17185
#>
#> $model_results$estimate_matrices$omega
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]
#> [1,] 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923
#> [2,] 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923
#> [,8] [,9] [,10] [,11] [,12]
#> [1,] 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923
#> [2,] 0.1142923 0.1142923 0.1142923 0.1142923 0.1142923
#>
#> $model_results$estimate_matrices$pdet
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]
#> [1,] 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842
#> [2,] 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842
#> [,8] [,9] [,10] [,11] [,12]
#> [1,] 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842
#> [2,] 0.9538842 0.9538842 0.9538842 0.9538842 0.9538842
#>
#>
#>
#> $model_data
#> $model_data$K
#> [1] 400
#>
#> $model_data$nit
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
#> [1,] 134 61 93 84 92 65 105 86 69 84 87 68
#> [2,] 148 93 90 103 94 58 65 65 53 63 67 63
#>
#> $model_data$l_s_c
#> NULL
#>
#> $model_data$g_s_c
#> NULL
#>
#> $model_data$g_t_c
#> NULL
#>
#> $model_data$o_s_c
#> NULL
#>
#> $model_data$o_t_c
#> NULL
#>
#> $model_data$p_s_c
#> NULL
#>
#> $model_data$p_t_c
#> NULL
#>
#> $model_data$SMALL_a_CORRECTION
#> [1] FALSE
Multi-threading is implemented using R packages doParallel and foreach. Note that multi-threading is used to split the computation of the transition probability matrix by rows. This may not be efficient on some architectures, and will not be efficient for small K: efficiency increases with increasing K. Alternative choices for multi-core processing are being considered.