Update nonparboot

10 files changed, 99 insertions, 234 deletions

diff --git a/R/boot.R b/R/boot.R
index 3da4fd0..202d8ef 100644
--- a/R/boot.R
+++ b/R/boot.R
@@ -153,11 +153,23 @@ setMethod("plot", signature(x="parboot", y="missing"),
 ##  they will be processed by vcov, confint, etc.
 
 setGeneric("nonparboot",
-    function(object, B = 0, ...) {standardGeneric("nonparboot")})
+  function(object, B = 0, keepOldSamples = TRUE, ...){
+  standardGeneric("nonparboot")
+})
+
+setMethod("nonparboot", "unmarkedFit", 
+  function(object, B = 0, keepOldSamples = TRUE, ...){
+  nonparboot_internal(object, B = B, keepOldSamples = keepOldSamples, ...)
+})
+
 
+setGeneric("nonparboot_internal", function(object, B = 0, keepOldSamples = TRUE, ...){
+  standardGeneric("nonparboot_internal")
+})
 
-setMethod("nonparboot", "unmarkedFit",
-          function(object, B = 0, keepOldSamples = TRUE, bsType, ...) {
+# "Site" method by default, other unmarkedFit types may set "both" below.
+setMethod("nonparboot_internal", "unmarkedFit",
+          function(object, B = 0, keepOldSamples = TRUE, bsType = "site", se = FALSE, ...) {
     bsType <- match.arg(bsType, c("site", "both"))
     if (identical(B, 0) && !is.null(object@bootstrapSamples)) {
         return(object)
@@ -185,7 +197,7 @@ setMethod("nonparboot", "unmarkedFit",
                           function(obs) sample(obs, replace = TRUE))
             data.b <- data.b[obs]
         }
-        fm <- update(object, data = data.b, se = FALSE)
+        fm <- update(object, data = data.b, se = se)
         return(fm)
     }
     if (!keepOldSamples) {
@@ -207,7 +219,7 @@ setMethod("nonparboot", "unmarkedFit",
 })
 
 
-setMethod("nonparboot", "unmarkedFitOccu",
+setMethod("nonparboot_internal", "unmarkedFitOccu",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
     callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
@@ -215,31 +227,7 @@ setMethod("nonparboot", "unmarkedFitOccu",
 })
 
 
-setMethod("nonparboot", "unmarkedFitPCount",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-setMethod("nonparboot", "unmarkedFitMPois",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-setMethod("nonparboot", "unmarkedFitDS",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-setMethod("nonparboot", "unmarkedFitOccuRN",
+setMethod("nonparboot_internal", "unmarkedFitOccuRN",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
     callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
@@ -247,36 +235,7 @@ setMethod("nonparboot", "unmarkedFitOccuRN",
 })
 
 
-
-setMethod("nonparboot", "unmarkedFitGMM",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-setMethod("nonparboot", "unmarkedFitGDS",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-
-
-setMethod("nonparboot", "unmarkedFitDailMadsen",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-
-
-
-setMethod("nonparboot", "unmarkedFitColExt",
+setMethod("nonparboot_internal", "unmarkedFitColExt",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
 #    browser()
@@ -339,63 +298,15 @@ setMethod("nonparboot", "unmarkedFitColExt",
 })
 
 
-setMethod("nonparboot", "unmarkedFitOccuPEN",
+setMethod("nonparboot_internal", "unmarkedFitOccuPEN",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
-#    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-#                   bsType="site")
-    bsType <- "site"
-    if (identical(B, 0) && !is.null(object@bootstrapSamples)) {
-        return(object)
-    }
-    if (B <= 0 && is.null(object@bootstrapSamples)) {
-        stop("B must be greater than 0 when fit has no bootstrap samples.")
-    }
-    data <- object@data
-    formula <- object@formula
-    designMats <- getDesign(data, formula) # bootstrap after removing sites
-    removed.sites <- designMats$removed.sites
-    if(length(removed.sites)>0)
-        data <- data[-removed.sites,]
-    y <- getY(data)
-    colnames(y) <- NULL
-    data@y <- y
-    M <- numSites(data)
-    boot.iter <- function() {
-        sites <- sort(sample(1:M, M, replace = TRUE))
-        data.b <- data[sites,]
-        y <- getY(data.b)
-        if (bsType == "both") {
-            obs.per.site <- lapply(1:nrow(y), function(i) which(!is.na(y[i,])))
-            obs <- lapply(obs.per.site,
-                          function(obs) sample(obs, replace = TRUE))
-            data.b <- data.b[obs]
-        }
-        fm <- update(object, data = data.b)
-        return(fm)
-    }
-    if (!keepOldSamples) {
-        object@bootstrapSamples <- NULL
-    }
-    object@bootstrapSamples <- c(object@bootstrapSamples,
-                                 replicate(B, boot.iter(),
-                                           simplify = FALSE))
-    coefs <- t(sapply(object@bootstrapSamples,
-                      function(x) coef(x)))
-    v <- cov(coefs)
-    object@covMatBS <- v
-    inds <- .estimateInds(object)
-    for (est in names(inds)) {
-        v.est <- v[inds[[est]], inds[[est]], drop = FALSE]
-        object@estimates@estimates[[est]]@covMatBS <- v.est
-    }
-    object
-
-
+    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
+                   bsType="site", se = TRUE)
 })
 
 
-setMethod("nonparboot", "unmarkedFitOccuPEN_CV",
+setMethod("nonparboot_internal", "unmarkedFitOccuPEN_CV",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
 #    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
@@ -456,16 +367,7 @@ setMethod("nonparboot", "unmarkedFitOccuPEN_CV",
 })
 
 
-setMethod("nonparboot", "unmarkedFitOccuTTD",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-
-setMethod("nonparboot", "unmarkedFitOccuMulti",
+setMethod("nonparboot_internal", "unmarkedFitOccuMulti",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
     bsType <- "site"
@@ -507,23 +409,6 @@ setMethod("nonparboot", "unmarkedFitOccuMulti",
     object
 })
 
-setMethod("nonparboot", "unmarkedFitNmixTTD",
-    function(object, B = 0, keepOldSamples = TRUE, ...)
-{
-    callNextMethod(object, B=B, keepOldSamples=keepOldSamples,
-                   bsType="site")
-})
-
-
-
-
-
-
-
-
-
-
-
 
 # ----------------------- Helper functions -------------------------------
 
diff --git a/R/gdistremoval.R b/R/gdistremoval.R
index 4610784..29c685f 100644
--- a/R/gdistremoval.R
+++ b/R/gdistremoval.R
@@ -827,7 +827,7 @@ setMethod("SSE", "unmarkedFitGDR", function(fit, ...){
     return(c(SSE = sum(r)))
 })
 
-setMethod("nonparboot", "unmarkedFitGDR",
+setMethod("nonparboot_internal", "unmarkedFitGDR",
     function(object, B = 0, keepOldSamples = TRUE, ...)
 {
    stop("Not currently supported for unmarkedFitGDR", call.=FALSE)
diff --git a/R/occuCOP.R b/R/occuCOP.R
index 672efd5..29ecc89 100644
--- a/R/occuCOP.R
+++ b/R/occuCOP.R
@@ -494,7 +494,7 @@ setMethod("simulate", "unmarkedFitOccuCOP",
 
 
 ## nonparboot ----
-setMethod("nonparboot", "unmarkedFitOccuCOP",
+setMethod("nonparboot_internal", "unmarkedFitOccuCOP",
   function(object, B = 0, keepOldSamples = TRUE, ...) {
   stop("Not currently supported for unmarkedFitOccuCOP", call.=FALSE)
 })
diff --git a/R/occuPEN.R b/R/occuPEN.R
index 55c7563..4fc0d50 100644
--- a/R/occuPEN.R
+++ b/R/occuPEN.R
@@ -168,9 +168,9 @@ occuPEN_CV <- function(formula, data, knownOcc = numeric(0), starts,
 
 occuPEN <- function(formula, data, knownOcc = numeric(0), starts,
                  method = "BFGS", engine = c("C", "R"),
-#		 se = TRUE,
 		 lambda = 0,
 		 pen.type = c("Bayes","Ridge","MPLE"),
+     se = FALSE,
 		 ...)
 {
 
diff --git a/man/nonparboot-methods.Rd b/man/nonparboot-methods.Rd
deleted file mode 100644
index 40dd72f..0000000
--- a/man/nonparboot-methods.Rd
+++ /dev/null
@@ -1,91 +0,0 @@
-\name{nonparboot-methods}
-\docType{methods}
-\alias{nonparboot}
-\alias{nonparboot-methods}
-\alias{nonparboot,unmarkedFit-method}
-\alias{nonparboot,unmarkedFitColExt-method}
-\alias{nonparboot,unmarkedFitDS-method}
-\alias{nonparboot,unmarkedFitMPois-method}
-\alias{nonparboot,unmarkedFitOccu-method}
-\alias{nonparboot,unmarkedFitOccuPEN-method}
-\alias{nonparboot,unmarkedFitOccuPEN_CV-method}
-\alias{nonparboot,unmarkedFitOccuRN-method}
-\alias{nonparboot,unmarkedFitPCount-method}
-\alias{nonparboot,unmarkedFitGDS-method}
-\alias{nonparboot,unmarkedFitGMM-method}
-\alias{nonparboot,unmarkedFitOccuTTD-method}
-\alias{nonparboot,unmarkedFitOccuMulti-method}
-\alias{nonparboot,unmarkedFitNmixTTD-method}
-\alias{nonparboot,unmarkedFitGDR-method}
-\alias{nonparboot,unmarkedFitDailMadsen-method}
-\alias{nonparboot,unmarkedFitOccuCOP-method}
-
-\title{ Nonparametric bootstrapping in unmarked }
-\description{
-Call \code{nonparboot} on an unmarkedFit to obtain non-parametric
-bootstrap samples.  These can then be used by \code{vcov} in order to
-get bootstrap estimates of standard errors.
-}
-\section{Methods}{
-\describe{
-
-\item{\code{signature(object = "unmarkedFit")}}{ Obtain nonparametric
-  bootstrap samples for a general unmarkedFit.  }
-
-\item{\code{signature(object = "unmarkedFitColExt")}}{ Obtain nonparametric
-  bootstrap samples for colext fits.  }
-
-\item{\code{signature(object = "unmarkedFitDS")}}{ Obtain nonparametric
-  bootstrap samples for a distsamp fits. }
-
-\item{\code{signature(object = "unmarkedFitMPois")}}{ Obtain nonparametric
-  bootstrap samples for a distsamp fits. }
-
-\item{\code{signature(object = "unmarkedFitOccu")}}{  Obtain nonparametric
-  bootstrap samples for a occu fits. }
-
-\item{\code{signature(object = "unmarkedFitOccuPEN")}}{  Obtain nonparametric
-  bootstrap samples for an occuPEN fit. }
-
-\item{\code{signature(object = "unmarkedFitOccuPEN_CV")}}{  Obtain nonparametric
-  bootstrap samples for occuPEN_CV fit. }
-
-\item{\code{signature(object = "unmarkedFitOccuRN")}}{  Obtain nonparametric
-  bootstrap samples for a occuRN fits. }
-
-\item{\code{signature(object = "unmarkedFitPCount")}}{  Obtain nonparametric
-  bootstrap samples for a pcount fits. }
-}}
-\details{
-
-  Calling \code{nonparboot} on an unmarkedFit returns the original
-  unmarkedFit, with the bootstrap samples added on.  Then subsequent
-  calls to \code{\link{vcov}} with the argument
-  \code{method="nonparboot"} will use these bootstrap samples.
-  Additionally, standard errors of derived estimates from either
-  \code{\link{linearComb}} or \code{\link{backTransform}} can be
-  instructed to use bootstrap samples by providing the argument
-  \code{method = "nonparboot"}.
-
-  For \code{\link{occu}} and \code{\link{occuRN}} both sites and
-  occassions are re-sampled. For all other fitting functions, only sites
-  are re-sampled.
-
-}
-\examples{
-data(ovendata)
-ovenFrame <- unmarkedFrameMPois(ovendata.list$data,
-siteCovs=as.data.frame(scale(ovendata.list$covariates[,-1])), type = "removal")
-(fm <- multinomPois(~ 1 ~ ufc + trba, ovenFrame))
-fm <- nonparboot(fm, B = 20) # should use larger B in real life.
-vcov(fm, method = "hessian")
-vcov(fm, method = "nonparboot")
-avg.abundance <- backTransform(linearComb(fm, type = "state", coefficients = c(1, 0, 0)))
-
-## Bootstrap sample information propagates through to derived quantities.
-vcov(avg.abundance, method = "hessian")
-vcov(avg.abundance, method = "nonparboot")
-SE(avg.abundance, method = "nonparboot")
-}
-\keyword{methods}
-
diff --git a/man/nonparboot.Rd b/man/nonparboot.Rd
new file mode 100644
index 0000000..7e1906d
--- /dev/null
+++ b/man/nonparboot.Rd
@@ -0,0 +1,60 @@
+\name{nonparboot}
+\alias{nonparboot}
+\alias{nonparboot-methods}
+\alias{nonparboot,unmarkedFit-method}
+
+\title{Nonparametric bootstrapping in unmarked}
+
+\description{
+Call \code{nonparboot} on an \code{unmarkedFit} to obtain non-parametric
+bootstrap samples.  These can then be used by \code{vcov} in order to
+get bootstrap estimates of standard errors. The function can be run
+repeatedly on an \code{unmarkedFit} object to add additional bootstrap samples.
+}
+
+\usage{
+\S4method{nonparboot}{unmarkedFit}(object, B = 0, keepOldSamples = TRUE, ...)
+}
+
+\arguments{
+  \item{object}{A \code{unmarkedFit} object}
+  \item{B}{Number of bootstrap samples to generate}
+  \item{keepOldSamples}{Logical. Should previously generated bootstrap samples
+        be kept?}
+  \item{...}{Other arguments, currently ignored} 
+}
+
+\value{An \code{unmarkedFit} object with added bootstrap samples.}
+
+\details{
+
+  Calling \code{nonparboot} on an unmarkedFit returns the original
+  unmarkedFit, with the bootstrap samples added on.  Then subsequent
+  calls to \code{\link{vcov}} with the argument
+  \code{method="nonparboot"} will use these bootstrap samples.
+  Additionally, standard errors of derived estimates from either
+  \code{\link{linearComb}} or \code{\link{backTransform}} can be
+  instructed to use bootstrap samples by providing the argument
+  \code{method = "nonparboot"}.
+
+  For \code{\link{occu}} and \code{\link{occuRN}} both sites and
+  occassions are re-sampled. For all other fitting functions, only sites
+  are re-sampled.
+
+}
+
+\examples{
+data(ovendata)
+ovenFrame <- unmarkedFrameMPois(ovendata.list$data,
+siteCovs=as.data.frame(scale(ovendata.list$covariates[,-1])), type = "removal")
+(fm <- multinomPois(~ 1 ~ ufc + trba, ovenFrame))
+fm <- nonparboot(fm, B = 20) # should use larger B in real life.
+vcov(fm, method = "hessian")
+vcov(fm, method = "nonparboot")
+avg.abundance <- backTransform(linearComb(fm, type = "state", coefficients = c(1, 0, 0)))
+
+## Bootstrap sample information propagates through to derived quantities.
+vcov(avg.abundance, method = "hessian")
+vcov(avg.abundance, method = "nonparboot")
+SE(avg.abundance, method = "nonparboot")
+}
diff --git a/man/occuPEN.Rd b/man/occuPEN.Rd
index 6b12d3e..7e14c20 100644
--- a/man/occuPEN.Rd
+++ b/man/occuPEN.Rd
@@ -5,7 +5,7 @@
 \title{Fit the MacKenzie et al. (2002) Occupancy Model with the penalized likelihood methods of Hutchinson et al. (2015)}
 
 \usage{occuPEN(formula, data, knownOcc=numeric(0), starts, method="BFGS",
-    engine=c("C", "R"), lambda=0, pen.type = c("Bayes","Ridge","MPLE"), ...)}
+    engine=c("C", "R"), lambda=0, pen.type = c("Bayes","Ridge","MPLE"), se = FALSE, ...)}
 
 \arguments{
     \item{formula}{Double right-hand side formula describing covariates of
@@ -20,6 +20,8 @@
       code during the optimization.}
     \item{lambda}{Penalty weight parameter.}
     \item{pen.type}{Which form of penalty to use.}
+    \item{se}{Nonfunctional; always FALSE; for compatibility with
+              other parts of \code{unmarked} }
     \item{\dots}{Additional arguments to optim, such as lower and upper
       bounds}
   }
diff --git a/tests/testthat/test_nonparboot.R b/tests/testthat/test_nonparboot.R
index d8fe127..1d4f539 100644
--- a/tests/testthat/test_nonparboot.R
+++ b/tests/testthat/test_nonparboot.R
@@ -63,6 +63,7 @@ test_that("nonparboot works with gmultmix", {
     fm1 <- gmultmix(~1, ~1, ~1, umf1, K=10)
     fm1 <- nonparboot(fm1, B=2)
     expect_equal(length(fm1@bootstrapSamples), 2)
+    expect_equal(vcov(fm1)[1,1], 0.1165165, tol=1e-7)
 
     umf2 <- unmarkedFrameGMM(y=y.ijt, siteCovs=data.frame(sc=sc1),
                              obsCovs=list(oc=oc1),
diff --git a/tests/testthat/test_occuMulti.R b/tests/testthat/test_occuMulti.R
index 6ce7921..85ff5d5 100644
--- a/tests/testthat/test_occuMulti.R
+++ b/tests/testthat/test_occuMulti.R
@@ -116,6 +116,11 @@ test_that("occuMulti can fit models with covariates",{
   detformulas <- c('~occ_cov1','~det_cov2')
   fm <- occuMulti(detformulas, stateformulas, data = umf, se=FALSE)
   expect_equivalent(coef(fm,'det')[2],3.355328e-05, tol=1e-4)
+
+  # nonparboot
+  set.seed(123)
+  npb <- nonparboot(fm, B=2)
+  expect_equal(vcov(npb, method="nonparboot")[1,1], 0.07733412, tol=1e-7)
 })
 
 test_that("occuMulti can handle NAs",{
diff --git a/tests/testthat/test_occuPEN.R b/tests/testthat/test_occuPEN.R
index 09e1d6c..2017559 100644
--- a/tests/testthat/test_occuPEN.R
+++ b/tests/testthat/test_occuPEN.R
@@ -120,10 +120,13 @@ test_that("occuPEN can fit models with covariates",{
   expect_error(fm <- occuPEN_CV(~ o1 + o2 ~ x, data = umf,foldAssignments=c(1,2,3,4,5),k=6))
 
   # nonparboot
+  set.seed(123)
   nbp <- nonparboot(fm, B=2)
   expect_is(nbp@covMatBS, "matrix")
+  expect_equal(nbp@covMatBS[1,1], 0.1691121, tol=1e-7)
   nbp_cv <- nonparboot(fmCV, B=2)
   expect_is(nbp_cv@covMatBS, "matrix")
+  expect_equal(nbp_cv@covMatBS[1,1], 0.5086074, tol=1e-7)
 
 })