summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--.gitignore1
-rw-r--r--src/plos-one.csl16
-rw-r--r--src/references.bib188
-rw-r--r--src/ubms-spatial.Rmd347
-rw-r--r--src/ubms-vignette.Rmd10
5 files changed, 557 insertions, 5 deletions
diff --git a/.gitignore b/.gitignore
index 43c785d..70e754d 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,3 +1,4 @@
*.html
*.xml
data/*
+*.Rds
diff --git a/src/plos-one.csl b/src/plos-one.csl
new file mode 100644
index 0000000..190a801
--- /dev/null
+++ b/src/plos-one.csl
@@ -0,0 +1,16 @@
+<?xml version="1.0" encoding="utf-8"?>
+<style xmlns="http://purl.org/net/xbiblio/csl" version="1.0" default-locale="en-US">
+ <!-- Generated with https://github.com/citation-style-language/utilities/tree/master/generate_dependent_styles/data/plos -->
+ <info>
+ <title>PLOS ONE</title>
+ <id>http://www.zotero.org/styles/plos-one</id>
+ <link href="http://www.zotero.org/styles/plos-one" rel="self"/>
+ <link href="http://www.zotero.org/styles/plos" rel="independent-parent"/>
+ <link href="http://www.plosone.org/static/guidelines#references" rel="documentation"/>
+ <category citation-format="numeric"/>
+ <category field="science"/>
+ <eissn>1932-6203</eissn>
+ <updated>2014-06-05T22:31:00+00:00</updated>
+ <rights license="http://creativecommons.org/licenses/by-sa/3.0/">This work is licensed under a Creative Commons Attribution-ShareAlike 3.0 License</rights>
+ </info>
+</style>
diff --git a/src/references.bib b/src/references.bib
new file mode 100644
index 0000000..d04362e
--- /dev/null
+++ b/src/references.bib
@@ -0,0 +1,188 @@
+@article{Bates_2015,
+ author = {Bates, Douglas and Mächler, Martin and Bolker, Ben and Walker, Steve},
+ title = {Fitting Linear Mixed-Effects Models Using lme4},
+ journal = {Journal of Statistical Software},
+ year = {2015},
+ number = {1},
+ volume = {67},
+ doi = {10.18637/jss.v067.i01}
+}
+
+@article{Carpenter_2017,
+ author = {Carpenter, Bob and Gelman, Andrew and Hoffman, Matthew D. and Lee, Daniel and Goodrich, Ben and Betancourt, Michael and Brubaker, Marcus and Guo, Jiqiang and Li, Peter and Riddell, Allen},
+ title = {Stan: A Probabilistic Programming Language},
+ journal = {Journal of Statistical Software},
+ publisher = {Foundation for Open Access Statistic},
+ year = {2017},
+ number = {1},
+ volume = {76},
+ url = {http://dx.doi.org/10.18637/jss.v076.i01},
+ doi = {10.18637/jss.v076.i01}
+}
+
+@article{Chandler2011,
+ author = {Chandler, Richard B. and Royle, J. Andrew and King, David I.},
+ title = {Inference about density and temporary emigration in unmarked populations},
+ journal = {Ecology},
+ year = {2011},
+ number = {7},
+ volume = {92},
+ pages = {1429--1435},
+ url = {http://dx.doi.org/10.1890/10-2433.1},
+ doi = {10.1890/10-2433.1},
+ file = {:Chandler2011b.pdf:PDF},
+ keyword = {Chestnut-sided Warbler, Dendroica pensylvanica, N-mixture model, USA, White Mountain National Forest, detection probability, hierarchical, models, population density, spot-mapping, temporary emigration, unmarked populations}
+}
+
+@article{Fiske_2011,
+ author = {Fiske, Ian and Chandler, Richard and others, },
+ title = {Unmarked: an {R} package for fitting hierarchical models of wildlife occurrence and abundance},
+ journal = {Journal of Statistical Software},
+ year = {2011},
+ number = {10},
+ volume = {43},
+ pages = {1--23},
+ doi = {10.18637/jss.v043.i10}
+}
+
+@book{Gelman2007,
+ author = {Gelman, Andrew and Hill, Jennifer},
+ title = {Data analysis using regression and multilevel/hierarchical models},
+ publisher = {Cambridge University Press},
+ year = {2007},
+ pages = {651},
+ url = {http://dx.doi.org/10.2277/0521867061, http://dx.doi.org/10.2277/0521867061},
+ doi = {10.2277/0521867061},
+ address = {New York, NY},
+ file = {:Gelman2007.pdf:PDF},
+ isbn = {0521867061},
+ pmid = {14341096}
+}
+
+@article{Hoffman_2014,
+ author = {Hoffman, Matthew D and Gelman, Andrew},
+ title = {The {N}o-{U}-{T}urn sampler: adaptively setting path lengths in {H}amiltonian {M}onte {C}arlo.},
+ journal = {Journal of Machine Learning Research},
+ year = {2014},
+ number = {1},
+ volume = {15},
+ pages = {1593--1623}
+}
+
+@article{MacKenzie_2002,
+ author = {MacKenzie, Darryl I. and Nichols, James D. and Lachman, Gideon B. and Droege, Sam and Royle, J. Andrew and Langtimm, Catherine A.},
+ title = {Estimating site occupancy rates when detection probabilities are less than one},
+ journal = {Ecology},
+ publisher = {Wiley},
+ year = {2002},
+ month = {aug},
+ number = {8},
+ volume = {83},
+ pages = {2248--2255},
+ url = {http://dx.doi.org/10.1890/0012-9658(2002)083[2248:esorwd]2.0.co;2},
+ doi = {10.1890/0012-9658(2002)083[2248:esorwd]2.0.co;2}
+}
+
+@article{MacKenzie_2004a,
+ author = {MacKenzie, Darryl I. and Bailey, Larissa L.},
+ title = {Assessing the fit of site-occupancy models},
+ journal = {Journal of Agricultural, Biological, and Environmental Statistics},
+ year = {2004},
+ number = {3},
+ volume = {9},
+ pages = {300--318},
+ doi = {10.1198/108571104x3361}
+}
+
+@article{Schmid_2004,
+ author = {Schmid, Hans and Zbinden, Niklaus and Keller, Verena},
+ title = {Überwachung der bestandsentwicklung häufiger brutvögel in der schweiz},
+ journal = {Swiss Ornithological Institute Sempach Switzerland},
+ year = {2004}
+}
+
+@article{Vehtari_2017,
+ author = {Vehtari, Aki and Gelman, Andrew and Gabry, Jonah},
+ title = {Practical {B}ayesian model evaluation using leave-one-out cross-validation and {WAIC}},
+ journal = {Statistics and Computing},
+ year = {2017},
+ number = {5},
+ volume = {27},
+ pages = {1413--1432},
+ doi = {10.1007/s11222-016-9696-4}
+}
+
+@article{Wright_2019,
+ author = {Wright, Wilson J. and Irvine, Kathryn M. and Higgs, Megan D.},
+ title = {Identifying occupancy model inadequacies: can residuals separately assess detection and presence?},
+ journal = {Ecology},
+ year = {2019},
+ pages = {e02703},
+ doi = {10.1002/ecy.2703}
+}
+
+@article{Yackulic_2020,
+ author = {Yackulic, Charles B. and Dodrill, Michael and Dzul, Maria and Sanderlin, Jamie S. and Reid, Janice A.},
+ title = {A need for speed in {B}ayesian population models: a practical guide to marginalizing and recovering discrete latent states},
+ journal = {Ecological Applications},
+ year = {2020},
+ number = {5},
+ volume = {30},
+ doi = {10.1002/eap.2112}
+}
+
+@article{Broms_2014,
+author = {Broms, Kristin M. and Johnson, Devin S. and Altwegg, Res and Conquest, Loveday L.},
+title = {Spatial occupancy models applied to atlas data show Southern Ground Hornbills strongly depend on protected areas},
+journal = {Ecological Applications},
+volume = {24},
+number = {2},
+pages = {363-374},
+doi = {https://doi.org/10.1890/12-2151.1},
+url = {https://esajournals.onlinelibrary.wiley.com/doi/abs/10.1890/12-2151.1},
+year = {2014}
+}
+
+
+@article{Clark_2019,
+author = {Clark, Allan E. and Altwegg, Res},
+title = {Efficient Bayesian analysis of occupancy models with logit link functions},
+journal = {Ecology and Evolution},
+volume = {9},
+number = {2},
+pages = {756-768},
+doi = {https://doi.org/10.1002/ece3.4850},
+url = {https://onlinelibrary.wiley.com/doi/abs/10.1002/ece3.4850},
+year = {2019}
+}
+
+
+@article{Hodges_2010,
+ title={Adding spatially-correlated errors can mess up the fixed effect you love},
+ author={Hodges, James S and Reich, Brian J},
+ journal={The American Statistician},
+ volume={64},
+ number={4},
+ pages={325--334},
+ year={2010},
+ publisher={Taylor \& Francis}
+}
+
+@article{Johnson_2013,
+author = {Johnson, Devin S. and Conn, Paul B. and Hooten, Mevin B. and Ray, Justina C. and Pond, Bruce A.},
+title = {Spatial occupancy models for large data sets},
+journal = {Ecology},
+volume = {94},
+number = {4},
+pages = {801-808},
+doi = {https://doi.org/10.1890/12-0564.1},
+url = {https://esajournals.onlinelibrary.wiley.com/doi/abs/10.1890/12-0564.1},
+year = {2013}
+}
+
+@book{Banerjee_2014,
+ title={Hierarchical modeling and analysis for spatial data},
+ author={Banerjee, Sudipto and Carlin, Bradley P and Gelfand, Alan E},
+ year={2014},
+ publisher={CRC press}
+}
diff --git a/src/ubms-spatial.Rmd b/src/ubms-spatial.Rmd
new file mode 100644
index 0000000..0d2d611
--- /dev/null
+++ b/src/ubms-spatial.Rmd
@@ -0,0 +1,347 @@
+---
+title: Spatial Models in ubms
+date: 2021-02-15
+bibliography: references.bib
+link-citations: yes
+csl: plos-one.csl
+output:
+ html_document:
+ css: style.css
+ highlight: pygments
+---
+
+# Introduction
+
+In many studies of animal and plant occurrence and abundance, sites may not be entirely independent.
+For example, sites may be close together, or site size may be smaller than a typical animal's home range.
+In these cases we expect spatial autocorrelation among nearby sites.
+Ideally, this autocorrelation would be accounted for in models.
+
+A solution to this problem is so-called spatial occupancy models, which include spatial random effects.
+A common approach to handling spatial autocorrelation is the use of an intrinsic conditional autoregressive (ICAR) random effect [@Banerjee_2014].
+The growing accessibility of Bayesian approaches has facilitated fitting these models.
+However, fitting ICAR models is computationally intensive and the ICAR random effect may be correlated with fixed covariates included in the model, making inference difficult [@Hodges_2010].
+To solve these issues, Johnson et al. [-@Johnson_2013] extended restricted spatial regression (RSR) to occupancy models.
+In an RSR model, the random effect is constructed in such a way that it is uncorrelated with the fixed covariates, and it can be significantly less computationally intensive than ICAR to estimate [@Broms_2014].
+Broms et al. [@Broms_2014] provide a nice overview of the concepts.
+Several recent studies have shown that RSR provides better and faster inference than ICAR [e.g. @Broms_2014; @Clark_2019].
+At least two specialized R packages ([`stocc`](https://cran.r-project.org/web/packages/stocc/index.html) and [`Rcppocc`](https://github.com/AllanClark/Rcppocc)) are available for fitting ICAR and RSR occupancy models.
+
+It is now possible to fit RSRs in Stan with `ubms` as well, allowing those familiar with `unmarked`/`ubms` syntax and workflows to easily fit spatial models.
+Currently only single-season occupancy models are supported, but support will soon be extended to other single-season model types.
+In addition to fitting RSRs, `ubms` also provides tools for applying and visualizing fitted spatial models.
+For now these features are available in a dev branch of `ubms` available on [Github](https://github.com/kenkellner/ubms/tree/RSR).
+
+# Example
+
+We will adapt the demo code and dataset available in R package `stocc` to fit a spatial occupancy model, which will also facilitate a comparison between `stocc` and `ubms` results.
+The `stocc` package was an extremely valuable guide when implementing RSRs in `ubms`, and contains many other useful features related to spatial occupancy modeling.
+
+```{r,eval=FALSE}
+library(ubms)
+library(stocc)
+```
+
+```{r, echo=FALSE, message=FALSE, warning=FALSE}
+suppressMessages(library(ubms))
+suppressMessages(library(stocc))
+```
+
+## Format the input data
+
+The sample dataset that comes with `stocc` is in long format.
+We need to convert it to the wide format that `unmarked` and `ubms` expect.
+
+
+```{r}
+# Load datasets
+data(habData) # site covariate data
+data(visitData) # observation cov data
+
+# Visits per site
+nobs <- table(visitData$site, visitData$obs)
+
+# Dimensions of y matrix
+J <- max(nobs) # max obs per site
+M <- nrow(habData) # no. of sites
+
+# Create blank matrices to hold wide-format versions
+y <- detCov1 <- detCov2 <- matrix(NA, nrow=M, ncol=J)
+
+# Iterate over sites
+for (i in 1:M){
+ # Data from site i
+ sub <- visitData[visitData$site==i,]
+ # Put this data into our wide matrices
+ subJ <- nrow(sub)
+ if(subJ == 0) next
+ y[i,1:subJ] <- sub$obs
+ detCov1[i,1:subJ] <- sub$detCov1
+ detCov2[i,1:subJ] <- as.character(sub$detCov2)
+}
+
+# Make site and obs cov data frames
+site_cov <- habData
+obs_cov <- data.frame(detCov1=as.vector(t(detCov1)), detCov2=as.vector(t(detCov2)))
+
+# Convert detCov2 back into a factor (as it was originally)
+obs_cov$detCov2 <- factor(obs_cov$detCov2, levels=c("0","1","2","3"))
+```
+
+In order to fit a spatial model, we need coordinates included in our site covariates.
+These should be in some projected coordinate system.
+In this dataset, coordinates are contained in columns `x` and `y`:
+
+```{r}
+head(site_cov)
+```
+
+Site and observation covariates are otherwise handled the same way as normal for `unmarked` and `ubms` models.
+Here's a quick visualization of how `habCov2` varies across the study area.
+
+```{r,fig.width=6}
+library(ggplot2)
+#Which sites were sampled at least once?
+sampled <- apply(y, 1, function(x) any(!is.na(x)))
+
+ggplot(data=site_cov, aes(x=x, y=y)) +
+ geom_tile(aes(fill=habCov2)) +
+ geom_point(data=site_cov[sampled,], size=0.5) +
+ scale_fill_gradientn(colors=heat.colors(10)) +
+ theme_bw(base_size=12) +
+ theme(panel.grid=element_blank())
+```
+
+The black dots represent locations that were sampled at least once.
+Note that not all grid cells with available site covariates were actually sampled; in fact most of them were not.
+For typical models, `ubms` would just ignore a site that was present in the dataset but never sampled (i.e., `y` all `NA`) since these sites cannot contribute any information.
+However for spatial models, generating predictive maps across areas that may not have been sampled is a typical goal.
+To facilitate this, spatial models in `ubms` do not discard sites which were never sampled but that have available site covariates.
+As you'd expect, you simply include these by including a row of all missing values in `y` and your observation covariates, as appropriate.
+For example, in this dataset, sites 1 and 2 were sampled but site 3 was not:
+
+```{r}
+y[1:3,]
+round(detCov1[1:3,],3)
+```
+
+The final step in data formatting is to construct the `unmarkedFrame`.
+
+```{r, warning=FALSE}
+umf <- unmarkedFrameOccu(y=y, siteCovs=site_cov, obsCovs=obs_cov)
+```
+
+## Choose RSR options
+
+Next, we need to decide on settings for the RSR.
+Specifically, we need to decide the distance threshold below which two sites will be considered neighbors and thus potentially correlated with each other.
+For example, if site A and site B are 0.5 distance units apart and we set the threshold at 1, A and B will be considered neighbors.
+If A and B are 1.5 distance units apart, they will not be neighbors.
+The distance units are defined by the units of the coordinates you provided.
+
+To visualize this, we can use the `RSR()` function, which we will use again later when fitting the model.
+The RSR() function requires as input at least one coordinate vector (typically you will have 2, e.g. `x` and `y`) and a value for the `threshold` argument.
+If we set the `plot_site` argument to an integer, `RSR` will return a map highlighting that site and its neighbors under the current threshold setting:
+
+```{r, fig.width=6}
+with(site_cov, RSR(x, y, threshold=1, plot_site=300))
+```
+
+With `threshold=1` and our coordinate grid, site 300 has just 4 neighbors.
+Now consider `threshold=10`:
+
+```{r, fig.width=6}
+with(site_cov, RSR(x, y, threshold=10, plot_site=300))
+```
+
+The threshold choice depends on your modeling goals and study system, but it is probably best to start with sites having a relatively small number of neighbors.
+We will set `threshold=1` for this example.
+
+A final argument available in `RSR()` is `moran_cut`.
+This option controls the number of eigenvectors that are used when calculating the spatial random effect (possible values are 1 to the number of sites).
+Generally, with fewer eigenvectors, the model will run faster and the result will be smoother.
+Previous studies have shown that inference is not particularly sensitive to this value, and recommended it be set to 10% of the number of sites [@Broms_2014; @Clark_2019].
+This is the default in `ubms` so we will leave it alone.
+
+## Fit the model
+
+We are finally ready to fit our spatial model.
+The formula syntax is identical to normal `ubms` and (and `unmarked`) models, except that we will include a call to `RSR()` with our desired options in our formula for `psi`.
+We also include two fixed detection and two fixed occupancy covariates.
+
+```{r}
+#Double formula: first part is for detection, second for occupancy
+form <- ~detCov1 + detCov2 ~habCov1 + habCov2 + RSR(x, y, threshold=1)
+```
+
+We will then fit the model in Stan with `stan_occu()`, first setting our desired number of parallel cores to use to 3.
+This will take Stan a few minutes with default MCMC settings.
+
+```{r, eval=FALSE}
+options(mc.cores=3)
+fit_ubms <- stan_occu(form, umf, chains=3)
+```
+
+```{r, message=FALSE, echo=FALSE, warning=FALSE}
+if(file.exists("fit_ubms.Rds")){
+ fit_ubms <- readRDS("fit_ubms.Rds")
+} else {
+ options(mc.cores=3)
+ fit_ubms <- stan_occu(form, umf, chains=3, refresh=0)
+ saveRDS(fit_ubms, "fit_ubms.Rds")
+}
+```
+
+You will likely get some warnings about effective sample size.
+The solution is to run the model for more iterations, but this is good enough for an example.
+
+## Examine results
+
+Look at the parameter estimates:
+
+```{r}
+fit_ubms
+```
+
+The precision term (tau, $\tau$) associated with the spatial random effect is the final row of the occupancy model estimates.
+Smaller values of $\tau$ imply greater variability in the spatial random effect.
+
+There also appears to be a strong effect of `habCov2` on occupancy.
+We can quickly visualize marginal effects of fixed covariates with `plot_marginal`.
+
+```{r, fig.width=5}
+plot_marginal(fit_ubms, "state")
+```
+
+You can extract predicted occupancy values for each site (including sites never actually sampled) using the `predict` function:
+
+```{r, message=FALSE}
+# "state" to get state/occupancy process, as opposed to "det"
+psi <- predict(fit_ubms, "state")
+head(psi)
+```
+
+## Plotting
+
+It is much more interesting to see how predicted occupancy varies across the study area.
+`ubms` has a built-in function for generating such a map, `plot_spatial`.
+
+```{r,message=FALSE, fig.width=6}
+plot_spatial(fit_ubms)
+```
+
+The plot also shows the location of sampled sites and whether the species was ever detected at each one.
+
+You can also look at the spatial distribution of random effects $\eta$:
+
+```{r,message=FALSE, fig.width=6}
+plot_spatial(fit_ubms, "eta")
+```
+
+## Model selection
+
+It is important to confirm that the spatial random effect is actually improving the predictive power of the model; if it isn't, we ought to remove it.
+To check this, first fit a similar model without the RSR:
+```{r, eval=FALSE}
+fit_nonspatial <- stan_occu(~detCov1 + detCov2 ~habCov1 + habCov2, umf, chains=3)
+```
+
+```{r, echo=FALSE, message=FALSE, warning=FALSE}
+if(file.exists("fit_nonspatial.Rds")){
+ fit_nonspatial <- readRDS("fit_nonspatial.Rds")
+} else {
+ options(mc.cores=3)
+ fit_nonspatial <- stan_occu(~detCov1 + detCov2 ~habCov1 + habCov2, umf,
+ chains=3, refresh=0)
+ saveRDS(fit_nonspatial, "fit_nonspatial.Rds")
+}
+```
+
+Combine the models to compare into a `fitList` (note: they must have the same number of chains and MCMC iterations):
+
+```{r}
+fl <- fitList(fit_ubms, fit_nonspatial)
+```
+
+Then, compare the spatial and non-spatial models with leave-one-out cross-validation (LOO)[@Vehtari_2017] using the `modSel()` function:
+
+```{r, warning=FALSE}
+round(modSel(fl),2)
+```
+
+A measure of predictive accuracy (`elpd`) is shown for each model, and the model with the highest predictive accuracy is listed first.
+The difference in `elpd` between the top and model and each subsequent model (`elpd_diff`), along with a measure of the uncertainty of that difference (`se_diff`), are also shown.
+If the magnitude of `elpd_diff` is large relative to the associated uncertainty, we can conclude that there is a difference in predictive power between the two models.
+In this case `elpd_diff` between the spatial (`fit_ubms`) and non-spatial (`fit_nonspatial`) models is several times the size of `se_diff`, so including the spatial random effect appears to improve model predictive accuracy.
+
+# Compare with stocc
+
+The implementation of spatial models in `ubms` owes much to the `stocc` package.
+The two packages give essentially identical results.
+To show this, we will fit the same model in `stocc` and compare the occupancy predictions with `ubms`.
+The code for fitting the model in `stocc` is taken directly from the package demo.
+
+```{r, eval=FALSE}
+names <- list(
+ visit=list(site="site",obs="obs"),
+ site=list(site="site", coords=c("x","y"))
+ )
+
+fit_stocc <- spatial.occupancy(
+ detection.model = ~ detCov1 + detCov2,
+ occupancy.model = ~ habCov1 + habCov2,
+ spatial.model = list(model="rsr", threshold=1, moran.cut=0.1*nrow(habData)),
+ so.data = make.so.data(visitData, habData, names),
+ prior = list(a.tau=0.5, b.tau=0.00005, Q.b=0.1, Q.g=0.1),
+ control = list(burnin=1000/5, iter=4000/5, thin=5)
+ )
+```
+
+```{r, echo=FALSE,message=FALSE, warning=FALSE}
+names <- list(
+ visit=list(site="site",obs="obs"),
+ site=list(site="site", coords=c("x","y"))
+ )
+
+if(file.exists("fit_stocc.Rds")){
+ fit_stocc <- readRDS("fit_stocc.Rds")
+} else {
+ # Takes a few minutes
+ fit_stocc <- spatial.occupancy(
+ detection.model = ~ detCov1 + detCov2,
+ occupancy.model = ~ habCov1 + habCov2,
+ spatial.model = list(model="rsr", threshold=1, moran.cut=0.1*nrow(habData)),
+ so.data = make.so.data(visitData, habData, names),
+ prior = list(a.tau=0.5, b.tau=0.00005, Q.b=0.1, Q.g=0.1),
+ control = list(burnin=1000/5, iter=4000/5, thin=5)
+ )
+ saveRDS(fit_stocc, "fit_stocc.Rds")
+}
+```
+
+Here's a similar plot of predicted occupancy across the study area:
+
+```{r, fig.width=6}
+# Predicted occupancy
+psi_stocc <- fit_stocc$occupancy.df$psi.est
+
+# Sites that were sampled at least once
+sampled <- apply(getY(umf), 1, function(x) any(!is.na(x)))
+plot_data <- cbind(habData, psi=psi_stocc)
+
+ggplot(data=plot_data, aes(x=x,y=y)) +
+ geom_point(aes(col=psi), size=5, pch=15) +
+ scale_color_gradientn(colors=terrain.colors(10)) +
+ geom_point(data=plot_data[sampled,],aes(x=x,y=y), size=0.5) +
+ labs(col="psi") +
+ theme_bw(base_size=12) +
+ theme(panel.grid=element_blank(),
+ strip.background=element_rect("white"))
+```
+
+As you can see, the predicted occupancy maps for `ubms` and `stocc` are almost identical.
+Note, however, that coefficient estimates may not be identical because `stocc` uses a probit link function, while `ubms` uses a logit link.
+
+# References
+
+<div id="refs"></div>
diff --git a/src/ubms-vignette.Rmd b/src/ubms-vignette.Rmd
index 8d328e3..870a6c7 100644
--- a/src/ubms-vignette.Rmd
+++ b/src/ubms-vignette.Rmd
@@ -22,7 +22,7 @@ knitr::opts_chunk$set(message=FALSE, warning=FALSE)
## What is ubms?
-`ubms` is an R package for fitting models of wildlife occurrence and abundance to with datasets where animals are not individually marked.
+`ubms` is an R package for fitting models of wildlife occurrence and abundance with datasets where animals are not individually marked.
It provides a nearly identical interface to the popular `unmarked` package [@Fiske_2011].
Instead of using maximum likelihood to fit models (as with `unmarked`),
models are fit in a Bayesian framework using [Stan](https://mc-stan.org/) [@Carpenter_2017].
@@ -150,7 +150,7 @@ options(mc.cores=3)
```{r, echo=FALSE}
library(ubms)
-(fit_stan <- stan_occu(~1~1, data=umf, chains=3, iter=500, refresh=0))
+(fit_stan <- stan_occu(~1~1, data=umf, chains=3, iter=1000, refresh=0))
```
### Compare results
@@ -232,13 +232,13 @@ fit_global <- stan_occu(~scale(date)~scale(forest)+scale(ele), data=umf,
```{r, echo=FALSE}
fit_null <- fit_stan
fit_sc <- stan_occu(~1~scale(forest)+scale(ele), data=umf,
- chains=3, iter=500, refresh=0)
-fit_oc <- stan_occu(~scale(date)~1, data=umf, chains=3, iter=500, refresh=0)
+ chains=3, iter=1000, refresh=0)
+fit_oc <- stan_occu(~scale(date)~1, data=umf, chains=3, iter=1000, refresh=0)
```
```{r, warning=TRUE, echo=FALSE}
fit_global <- stan_occu(~scale(date)~scale(forest)+scale(ele), data=umf,
- chains=3, iter=500, refresh=0)
+ chains=3, iter=1000, refresh=0)
```
The `fit_global` model gave us some warnings about the effective sample size (`n_eff`) along with a suggested solution.