diff options
author | Ken Kellner <ken@kenkellner.com> | 2022-07-22 17:09:21 -0400 |
---|---|---|
committer | Ken Kellner <ken@kenkellner.com> | 2022-07-22 17:09:21 -0400 |
commit | 288458c5e7e3fb5775c21d22ac8e4a1c63cbc206 (patch) | |
tree | f858b2da36824d2209f4c9303882f9d7379e737c | |
parent | e8c23369e3d2670e1b632abf818c572c69023071 (diff) |
Fix some bugs and add tests
-rw-r--r-- | DESCRIPTION | 4 | ||||
-rw-r--r-- | R/predict.R | 21 | ||||
-rw-r--r-- | R/simulate.R | 4 | ||||
-rw-r--r-- | tests/testthat/test_distsamp.R | 7 | ||||
-rw-r--r-- | tests/testthat/test_multinomPois.R | 6 | ||||
-rw-r--r-- | tests/testthat/test_powerAnalysis.R | 11 | ||||
-rw-r--r-- | tests/testthat/test_simulate.R | 9 |
7 files changed, 55 insertions, 7 deletions
diff --git a/DESCRIPTION b/DESCRIPTION index fe34bfc..972f834 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: unmarked -Version: 1.2.5.9002 -Date: 2022-06-17 +Version: 1.2.5.9003 +Date: 2022-07-22 Type: Package Title: Models for Data from Unmarked Animals Authors@R: c( diff --git a/R/predict.R b/R/predict.R index cd14637..68efbb7 100644 --- a/R/predict.R +++ b/R/predict.R @@ -132,9 +132,12 @@ setGeneric("get_formula", function(object, type, ...){ }) setMethod("get_formula", "unmarkedFit", function(object, type, ...){ - if(type == "state") return(as.formula(paste("~", object@formula[3], sep=""))) - if(type == "det") return(as.formula(object@formula[[2]])) - stop("Invalid type") + if(type == "state"){ + return(as.formula(paste("~", object@formula[3], sep=""))) + } else if(type == "det"){ + return(as.formula(object@formula[[2]])) + } + NULL }) # When newdata is data.frame/raster, get original dataset @@ -574,6 +577,12 @@ setMethod("get_orig_data", "unmarkedFitGDR", function(object, type, ...){ # bespoke predict method since it has numerious unusual options # and requires bootstrapping +# This method is used by simulate but not by predict +setMethod("get_formula", "unmarkedFitOccuMulti", function(object, type, ...){ + switch(type, state=object@stateformulas, + det=object@detformulas) +}) + setMethod("predict", "unmarkedFitOccuMulti", function(object, type, newdata, #backTransform = TRUE, na.rm = TRUE, @@ -748,6 +757,12 @@ setMethod("predict", "unmarkedFitOccuMulti", # bespoke predict method since it requires bootstrapping +# This method is used by simulate by not by predict +setMethod("get_formula", "unmarkedFitOccuMS", function(object, type, ...){ + switch(type, psi=object@psiformulas, phi=object@phiformulas, + det=object@detformulas) +}) + setMethod("predict", "unmarkedFitOccuMS", function(object, type, newdata, #backTransform = TRUE, na.rm = TRUE, diff --git a/R/simulate.R b/R/simulate.R index db10552..5df7247 100644 --- a/R/simulate.R +++ b/R/simulate.R @@ -109,9 +109,9 @@ generate_random_effects <- function(coefs, fit){ } if(!is.factor(lvldata)){ - stop("Random effect covariates must be specified as factors with guide argument", call.=FALSe) + stop("Random effect covariates must be specified as factors with guide argument", call.=FALSE) } - b <- rnorm(length(levels(lvldata)), 0, old_coefs[signame]) + b <- stats::rnorm(length(levels(lvldata)), 0, old_coefs[signame]) names(b) <- rep(paste0("b_",i), length(b)) new_coefs <- c(new_coefs, b) coefs[[i]] <- new_coefs diff --git a/tests/testthat/test_distsamp.R b/tests/testthat/test_distsamp.R index 9c42119..9d6fe04 100644 --- a/tests/testthat/test_distsamp.R +++ b/tests/testthat/test_distsamp.R @@ -265,6 +265,7 @@ test_that("getP works with distsamp",{ test_that("distsamp works with random effects",{ + set.seed(123) data(linetran) umf <- unmarkedFrameDS(y=as.matrix(linetran[,1:4]), siteCovs=linetran[,6:7], survey="line", tlength=linetran$Length, unitsIn='m', @@ -301,4 +302,10 @@ test_that("distsamp works with random effects",{ pr <- lapply(mods, function(x) predict(x, "state")) expect_true(all(sapply(pr, inherits, "data.frame"))) + # Make sure simulate accounts for random effects + s <- simulate(hn, nsim=30) + avg <- apply(sapply(s, function(x) apply(x,1,sum)),1, mean) + # average first count and predicted abundance should be highly correlated + pr <- predict(hn, 'state') + expect_true(cor(avg, pr$Predicted) > 0.7) }) diff --git a/tests/testthat/test_multinomPois.R b/tests/testthat/test_multinomPois.R index ee51335..6b4b693 100644 --- a/tests/testthat/test_multinomPois.R +++ b/tests/testthat/test_multinomPois.R @@ -237,6 +237,12 @@ test_that("multinomPois can fit models with random effects",{ expect_equivalent(dim(pr), c(100, 4)) expect_equivalent(dim(pr2), c(5,4)) + # Make sure simulate accounts for random effects + s <- simulate(fm, nsim=30) + avg <- apply(sapply(s, function(x) x[,1]),1, mean) + # average first count and predicted abundance should be highly correlated + expect_true(cor(avg, pr$Predicted) > 0.7) + umf2@y[1,1] <- NA umf2@y[2,] <- NA umf2@siteCovs$x1[3] <- NA diff --git a/tests/testthat/test_powerAnalysis.R b/tests/testthat/test_powerAnalysis.R index fbf75a4..d316a64 100644 --- a/tests/testthat/test_powerAnalysis.R +++ b/tests/testthat/test_powerAnalysis.R @@ -61,6 +61,17 @@ test_that("powerAnalysis method works",{ pl <- unmarkedPowerList(template_model, effect_sizes, design=scenarios, nsim=10) expect_is(pl, "unmarkedPowerList") + # With random effect + set.seed(123) + rguide <- list(group=factor(levels=letters[1:20])) + rform <- list(state=~x+(1|group), det=~1) + rcf <- list(state=c(intercept=0, x=0.5, group=0.7), det=c(intercept=0)) + umfr <- simulate("occu", formulas=rform, design=design, coefs=rcf, guide=rguide) + fm <- occu(~1~x+(1|group), umfr) + pa5 <- powerAnalysis(fm, rcf, nsim=10) + s <- summary(pa5) + expect_equal(nrow(s), 3) + expect_equal(s$Power[2], 1) }) }) diff --git a/tests/testthat/test_simulate.R b/tests/testthat/test_simulate.R index c1cda3e..04f1686 100644 --- a/tests/testthat/test_simulate.R +++ b/tests/testthat/test_simulate.R @@ -34,6 +34,15 @@ test_that("simulate can generate new datasets from scratch",{ expect_true(is.factor(umf2@siteCovs$landcover)) expect_equivalent(mean(umf2@siteCovs$elev), 2.01722, tol=1e-5) + # With random effect + set.seed(123) + rguide <- list(group=factor(levels=letters[1:20])) + rform <- list(state=~(1|group), det=~1) + rcf <- list(state=c(intercept=0, group=0.7), det=c(intercept=0)) + umfr <- simulate("occu", formulas=rform, design=design, coefs=rcf, guide=rguide) + fm <- occu(~1~(1|group), umfr) + expect_equal(sigma(fm)$sigma, 0.6903913, tol=1e-5) + # pcount set.seed(123) cf$alpha <- c(alpha=0.5) |