Fix problem with occuMS NA handling and a typo in docs

6 files changed, 39 insertions, 15 deletions

diff --git a/DESCRIPTION b/DESCRIPTION
index 972f834..751cc40 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: unmarked
-Version: 1.2.5.9003
-Date: 2022-07-22
+Version: 1.2.5.9004
+Date: 2022-07-27
 Type: Package
 Title: Models for Data from Unmarked Animals
 Authors@R: c(
diff --git a/R/ranef.R b/R/ranef.R
index 5547214..9dca184 100644
--- a/R/ranef.R
+++ b/R/ranef.R
@@ -121,6 +121,7 @@ setMethod("ranef", "unmarkedFitOccuMS", function(object, ...)
       g <- rep(1, S)
       p_raw <- sapply(p_all, function(x) x[i,])
       for (j in 1:nrow(p_raw)){
+        if(any(is.na(p_raw[j,])) | is.na(y[i,j])) next
         sdp <- matrix(0, nrow=S, ncol=S)
         sdp[guide] <- p_raw[j,]
         sdp[,1] <- 1 - rowSums(sdp)
@@ -142,6 +143,7 @@ setMethod("ranef", "unmarkedFitOccuMS", function(object, ...)
       p_raw <- sapply(p_all, function(x) x[i,])
       for (j in 1:nrow(p_raw)){
         probs <- p_raw[j,]
+        if(any(is.na(probs)) | is.na(y[i,j])) next
         sdp <- matrix(0, nrow=S, ncol=S)
         sdp[1,1] <- 1
         sdp[2,1:2] <- c(1-probs[1], probs[1])
diff --git a/R/unmarkedFit.R b/R/unmarkedFit.R
index 6557b9d..a672c6f 100644
--- a/R/unmarkedFit.R
+++ b/R/unmarkedFit.R
@@ -2552,14 +2552,15 @@ setMethod("simulate", "unmarkedFitOccuMS",
   for (n in 1:N){
     yindex <- 1
     for (t in 1:T){
-      if (z[n,t] == 0) {
-        yindex <- yindex + J
-        next
-      }
       for (j in 1:J){
-
         if(prm == "multinomial"){
           probs_raw <- sapply(p, function(x) x[n,yindex])
+          # Make sure output is NA if probs have NA
+          if(any(is.na(probs_raw))){
+            y[n,yindex] <- NA
+            yindex <- yindex + 1
+            next
+          }
 
           sdp <- matrix(0, nrow=S, ncol=S)
           sdp[guide] <- probs_raw
@@ -2571,13 +2572,22 @@ setMethod("simulate", "unmarkedFitOccuMS",
           p11 <- p[[1]][n,yindex]
           p12 <- p[[2]][n,yindex]
           p22 <- p[[3]][n,yindex]
+          # Trap NAs in probability of detection
+          if(any(is.na(c(p11, p12, p22)))){
+            y[n,yindex] <- NA
+            next
+          }
           probs <- switch(z[n,t]+1,
                           c(1,0,0),
                           c(1-p11,p11,0),
                           c(1-p12,p12*(1-p22),p12*p22))
         }
-
-        y[n,yindex] <- sample(0:(S-1), 1, prob=probs)
+        # this NA trap probably isn't necessary but leaving it in just in case
+        if(all(!is.na(probs))){
+          y[n,yindex] <- sample(0:(S-1), 1, prob=probs)
+        } else {
+          y[n,yindex] <- NA
+        }
         yindex <- yindex + 1
       }
     }
diff --git a/man/formatDistData.Rd b/man/formatDistData.Rd
index e4b1b66..f5aaa83 100644
--- a/man/formatDistData.Rd
+++ b/man/formatDistData.Rd
@@ -19,7 +19,7 @@ transect names.}
   than once, this can be used to format data for \code{gdistsamp}. It is
   the name of the column in distData that contains the occasion
   numbers. The occasion column should be a factor.}
-\item{effortMatrix}{optional matrix of 1 and 0s that is M * J in size and will allow for the insertion of NAs where the matrix = 0, indicating that a survey was not completed. When not supplied a matrix of all 1s is created since it is assumed all surveys were completed.}
+\item{effortMatrix}{optional matrix of 1 and 0s that is M * T in size and will allow for the insertion of NAs where the matrix = 0, indicating that a survey was not completed. When not supplied a matrix of all 1s is created since it is assumed all surveys were completed.}
 }
 \details{This function creates a site (M) by distance interval (J) response
 matrix from a data.frame containing the detection distances for each
diff --git a/tests/testthat/test_occuMS.R b/tests/testthat/test_occuMS.R
index 17ca610..aa0d024 100644
--- a/tests/testthat/test_occuMS.R
+++ b/tests/testthat/test_occuMS.R
@@ -213,7 +213,7 @@ test_that("occuMS can fit the multinomial model",{
   expect_equivalent(length(sim),3)
   expect_true(all(unlist(sim)%in%c(0:2)))
   expect_equivalent(mean(fit_C@data@y),0.268)
-  expect_equivalent(sapply(sim,mean),c(0.244,0.280,0.288))
+  expect_equivalent(sapply(sim,mean),c(0.232,0.252,0.276))
 
   #check fitted
   set.seed(123)
@@ -326,7 +326,7 @@ test_that("occuMS can fit the conditional binomial model",{
   expect_equivalent(length(sim),3)
   expect_true(all(unlist(sim)%in%c(0:2)))
   expect_equivalent(mean(fit_C@data@y),0.2)
-  expect_equivalent(sapply(sim,mean),c(0.200,0.156,0.128))
+  expect_equivalent(sapply(sim,mean),c(0.172,0.196,0.184))
 })
 
 test_that("occuMS handles NAs properly",{
@@ -361,10 +361,22 @@ test_that("occuMS handles NAs properly",{
 
   yna <- y
   yna[1,1] <- NA
+  obs_covs[1,1] <- NA
   umf <- unmarkedFrameOccuMS(y=yna,siteCovs=site_covs,obsCovs=obs_covs)
   fit <- occuMS(rep('~1',3),rep('~1',2),data=umf,se=F)
   expect_equivalent(fit@AIC,53.06711,tol=1e-4)
 
+  # Check simulate and ranef methods
+  fit <- occuMS(rep('~V1',3),rep('~1',2),data=umf,se=F)
+  s <- simulate(fit, nsim=3)
+  expect_equal(sum(is.na(unlist(s))), 3)
+  r <- ranef(fit)
+  expect_true(!any(is.na(r@post)))
+
+  fit_cb <- occuMS(rep('~V1',3),rep('~1',2),data=umf,se=F, parameterization='condbinom')
+  s <- simulate(fit_cb, nsim=3)
+  expect_equal(sum(is.na(unlist(s))), 3)
+
   yna <- y
   yna[1,] <- NA
   sc_na <- site_covs
@@ -514,7 +526,7 @@ expect_equivalent(length(coef(fit_new)),14)
 
 set.seed(123)
 fit_sim <- simulate(fitC,nsim=2)
-expect_equivalent(fit_sim[[1]][2,],c(0,2,1,0,0,2))
+expect_equivalent(fit_sim[[1]][2,],c(0,0,0,0,0,0))
 
 nul <- capture.output(pr_phi <- predict(fitC,'phi'))
 pr_phi <- sapply(pr_phi, function(x) x$Predicted[1])
@@ -642,7 +654,7 @@ expect_equivalent(fit_cbC@AIC,820.0645,tol=1e-4)
 
 set.seed(123)
 fit_sim <- simulate(fit_cbC,nsim=1)
-expect_equivalent(fit_sim[[1]][1,],c(0,0,0,2,1,0))
+expect_equivalent(fit_sim[[1]][1,],c(0,0,0,0,2,1))
 
 nul <- capture.output(pr_phi <- predict(fit_cbC,'phi'))
 pr_phi <- sapply(pr_phi, function(x) x$Predicted[1])
diff --git a/tests/testthat/test_simulate.R b/tests/testthat/test_simulate.R
index 04f1686..be8a356 100644
--- a/tests/testthat/test_simulate.R
+++ b/tests/testthat/test_simulate.R
@@ -182,7 +182,7 @@ test_that("simulate can generate new datasets from scratch",{
   cf <- list(state=bstate, det=bdet)
   expect_warning(umf15 <- simulate("occuMS", formulas=forms, coefs=cf, design=list(M=500, J=5, T=1)))
   fm <- occuMS(forms$det, forms$state, data=umf15, parameterization="multinomial")
-  expect_equivalent(coef(fm, 'state'), c(-0.437,0.767,-0.671,-0.595), tol=1e-3)
+  expect_equivalent(coef(fm, 'state'), c(-0.657,1.033,-0.633,-0.582), tol=1e-3)
 
   # gdistremoval
   set.seed(123)