--- title: "Data-augmented models in flocker" author: "Jacob Socolar" date: "2023-10-20" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Data-augmented models in flocker} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- When modeling multiple ecologically comparable species simultaneously, occupancy models often assume conditional exchangeability across species and fit species-specific terms as random effects. Data-augmented multi-species models leverage this exchangeability assumption to estimate the number and prevalence of never-detected species within the study region. To do so, the dataset is augmented with a large number of pseudospecies with all-zero detection histories, each species (both the observed species and the augmented pseudospecies) is ascribed a common parameter $\omega$ giving the Bernoulli probability that a given (pseudo)species occurs in the study area, and random-effects exchangeability assumptions are assumed to hold for all species-specific modeled terms (i.e. intercepts and slopes for occupancy and detection). For a data-augmented multi-species model, we marginalize over the occupancy status of a closure-unit as for a single-season model yielding the unit-wise likelihood $\mathcal{L}_i$, and we additionally marginalize over the availability of each (pseudo)species, yielding the species-wise likelihood \[ \mathcal{N}_s = \begin{cases} B(1 | \omega)\prod\limits_{i \textrm{ in } I_s}{\mathcal{L}_i} & \text{if $r_s = 1$} \\ B(0 | \omega) + B(1 | \omega)\prod\limits_{i \textrm{ in } I_s}{\mathcal{L}_i}& \text{if $r_s = 0$} \end{cases} \] where $s$ indexes the species, $\omega$ is the fitted availability probability of a species, $I_s$ is the set of all closure-unit indices $i$ pertaining to species $s$, $r_s$ is an indicator that takes the value $1$ if there is at least one positive detection of the (pseudo)species in the entire dataset and $0$ otherwise. Fitting the data-augmented model in `flocker` requires passing the observed data as a three-dimensional array with sites along the first dimension, visits along the second, and species along the third. Additionally, we must supply the `n_aug` argument to `make_flocker_data()`, specifying how many all-zero pseudospecies to augment the data with. ``` #> Formatting data for a single-season multispecies occupancy model with data augmentation for never-observed species. For details, see make_flocker_data_augmented. All warnings and error messages should be interpreted in the context of make_flocker_data_augmented #> formatting rep indices #> | | | 0% | | | 1% | |= | 1% | |= | 2% | |== | 2% | |== | 3% | |=== | 3% | |=== | 4% | |=== | 5% | |==== | 5% | |==== | 6% | |===== | 6% | |===== | 7% | |====== | 7% | |====== | 8% | |====== | 9% | |======= | 9% | |======= | 10% | |======== | 10% | |======== | 11% | |========= | 11% | |========= | 12% | |========== | 12% | |========== | 13% | |========== | 14% | |=========== | 14% | |=========== | 15% | |============ | 15% | |============ | 16% | |============= | 16% | |============= | 17% | |============= | 18% | |============== | 18% | |============== | 19% | |=============== | 19% | |=============== | 20% | |================ | 20% | |================ | 21% | |================ | 22% | |================= | 22% | |================= | 23% | |================== | 23% | |================== | 24% | |=================== | 24% | |=================== | 25% | |=================== | 26% | |==================== | 26% | |==================== | 27% | |===================== | 27% | |===================== | 28% | |====================== | 28% | |====================== | 29% | |====================== | 30% | |======================= | 30% | |======================= | 31% | |======================== | 31% | |======================== | 32% | |========================= | 32% | |========================= | 33% | |========================= | 34% | |========================== | 34% | |========================== | 35% | |=========================== | 35% | |=========================== | 36% | |============================ | 36% | |============================ | 37% | |============================ | 38% | |============================= | 38% | |============================= | 39% | |============================== | 39% | |============================== | 40% | |=============================== | 40% | |=============================== | 41% | |================================ | 41% | |================================ | 42% | |================================ | 43% | |================================= | 43% | |================================= | 44% | |================================== | 44% | |================================== | 45% | |=================================== | 45% | |=================================== | 46% | |=================================== | 47% | |==================================== | 47% | |==================================== | 48% | |===================================== | 48% | |===================================== | 49% | |====================================== | 49% | |====================================== | 50% | |====================================== | 51% | |======================================= | 51% | |======================================= | 52% | |======================================== | 52% | |======================================== | 53% | |========================================= | 53% | |========================================= | 54% | |========================================= | 55% | |========================================== | 55% | |========================================== | 56% | |=========================================== | 56% | |=========================================== | 57% | |============================================ | 57% | |============================================ | 58% | |============================================ | 59% | |============================================= | 59% | |============================================= | 60% | |============================================== | 60% | |============================================== | 61% | |=============================================== | 61% | |=============================================== | 62% | |================================================ | 62% | |================================================ | 63% | |================================================ | 64% | |================================================= | 64% | |================================================= | 65% | |================================================== | 65% | |================================================== | 66% | |=================================================== | 66% | |=================================================== | 67% | |=================================================== | 68% | |==================================================== | 68% | |==================================================== | 69% | |===================================================== | 69% | |===================================================== | 70% | |====================================================== | 70% | |====================================================== | 71% | |====================================================== | 72% | |======================================================= | 72% | |======================================================= | 73% | |======================================================== | 73% | |======================================================== | 74% | |========================================================= | 74% | |========================================================= | 75% | |========================================================= | 76% | |========================================================== | 76% | |========================================================== | 77% | |=========================================================== | 77% | |=========================================================== | 78% | |============================================================ | 78% | |============================================================ | 79% | |============================================================ | 80% | |============================================================= | 80% | |============================================================= | 81% | |============================================================== | 81% | |============================================================== | 82% | |=============================================================== | 82% | |=============================================================== | 83% | |=============================================================== | 84% | |================================================================ | 84% | |================================================================ | 85% | |================================================================= | 85% | |================================================================= | 86% | |================================================================== | 86% | |================================================================== | 87% | |================================================================== | 88% | |=================================================================== | 88% | |=================================================================== | 89% | |==================================================================== | 89% | |==================================================================== | 90% | |===================================================================== | 90% | |===================================================================== | 91% | |====================================================================== | 91% | |====================================================================== | 92% | |====================================================================== | 93% | |======================================================================= | 93% | |======================================================================= | 94% | |======================================================================== | 94% | |======================================================================== | 95% | |========================================================================= | 95% | |========================================================================= | 96% | |========================================================================= | 97% | |========================================================================== | 97% | |========================================================================== | 98% | |=========================================================================== | 98% | |=========================================================================== | 99% | |============================================================================| 99% | |============================================================================| 100% #> Compiling Stan program... #> Start sampling #> Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable. #> Running the chains for more iterations may help. See #> https://mc-stan.org/misc/warnings.html#bulk-ess ``` Here, the random effect of species is specified using the special grouping keyword `ff_species` (names beginning with `ff_` are reserved in `flocker` and are not allowed as names for user-supplied covariates). ```r summary(fm) #> Family: occupancy_augmented #> Links: mu = identity; occ = identity; Omega = identity #> Formula: ff_y | vint(ff_n_unit, ff_n_rep, ff_Q, ff_n_sp, ff_superQ, ff_species, ff_rep_index1, ff_rep_index2, ff_rep_index3, ff_rep_index4) ~ uc1 + ec1 + (1 + uc1 + ec1 | ff_species) #> occ ~ (1 | ff_species) #> Omega ~ 1 #> Data: data (Number of observations: 25600) #> Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1; #> total post-warmup draws = 4000 #> #> Group-Level Effects: #> ~ff_species (Number of levels: 128) #> Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS #> sd(Intercept) 1.59 0.26 1.16 2.17 1.00 753 1572 #> sd(uc1) 0.25 0.11 0.06 0.48 1.00 1303 1041 #> sd(ec1) 0.39 0.09 0.23 0.60 1.00 2010 2878 #> sd(occ_Intercept) 1.59 0.32 1.09 2.36 1.01 881 1406 #> cor(Intercept,uc1) 0.74 0.22 0.14 0.97 1.00 2155 1983 #> cor(Intercept,ec1) 0.57 0.19 0.11 0.85 1.00 2915 2617 #> cor(uc1,ec1) 0.58 0.26 -0.05 0.94 1.01 913 1378 #> #> Population-Level Effects: #> Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS #> Intercept 0.47 0.32 -0.18 1.09 1.01 368 818 #> occ_Intercept 0.21 0.36 -0.52 0.87 1.01 454 983 #> Omega_Intercept -1.26 0.22 -1.70 -0.85 1.00 6322 2762 #> uc1 0.01 0.08 -0.14 0.16 1.00 1061 1985 #> ec1 0.07 0.09 -0.11 0.25 1.00 1479 2479 #> #> Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS #> and Tail_ESS are effective sample size measures, and Rhat is the potential #> scale reduction factor on split chains (at convergence, Rhat = 1). ``` `flocker` enables users to fit data-augmented models using arbitrary `brms` formulas for the occupancy and detection components. However, we caution that continuous covariates in the occupancy sub-model can lead to pitfalls in interpretation. A seemingly straightforward application, for example, might be to ask how many species are present along an elevational gradient, fitting species-specific quadratic elevation-occupancy relationships. In our experience, a data-augmented model that includes quadratic elevation terms will place an arbitrarily large number of pseudo-species along the gradient, but will do so by placing pseudospecies' estimated elevational ranges entirely outside the range covered by the sampling effort ([Socolar et al 2022](https://onlinelibrary.wiley.com/doi/10.1002/ece3.9328)). The model is in effect trying to estimate how many species occur in a landscape with elevations ranging from negative to positive infinity. This extrapolation is unprincipled, most obviously because it does not account for hard limits imposed by the physical termina of the gradient (valley floor and mountain peak). Thus, although `flocker` provides functionality to fit continuous covariates in the occupancy term, we recommend extreme caution in interpreting patterns estimated for never-observed species.