testthat::test_that("qvcalc.default works for glm", {
## Overdispersed Poisson loglinear model for ship damage data
## from McCullagh and Nelder (1989), Sec 6.3.2
testthat::skip_if_not_installed("MASS")
require(MASS)
data(ships)
ships$year <- as.factor(ships$year)
ships$period <- as.factor(ships$period)
shipmodel <- glm(formula = incidents ~ type + year + period,
family = quasipoisson,
data = ships,
subset = (service > 0),
offset = log(service))
shiptype.qv <- qvcalc(shipmodel, "type")
testthat::expect_named(
shiptype.qv,
c("covmat", "qvframe", "dispersion", "relerrs", "factorname",
"coef.indices", "modelcall")
)
testthat::expect_equal(
shiptype.qv$qvframe$quasiSE,
c(0.200965853501369, 0.112705161723464,
0.375250492322657, 0.323893536427565,
0.232171862907782)
)
testthat::expect_error({
qvcalc(shipmodel, "period")
}, regexp = "factors with 3 or more levels")
testthat::expect_error({
qvcalc(shipmodel)
}, regexp = "both NULL")
summary(shiptype.qv)
print(shiptype.qv)
testthat::expect_silent({
plot(shiptype.qv)
})
# no intercept model
shipmodel2 <- glm(formula = incidents ~ type + year + period - 1,
family = quasipoisson,
data = ships,
subset = (service > 0),
offset = log(service))
shiptype.qv2 <- qvcalc(shipmodel2, "type")
testthat::expect_equal(
shiptype.qv$qvframe$quasiSE,
shiptype.qv2$qvframe$quasiSE
)
# no intercept model - using coef indices
shiptype.qv2 <- qvcalc(shipmodel2, coef.indices = 1:5)
testthat::expect_equal(
shiptype.qv$qvframe$quasiSE,
shiptype.qv2$qvframe$quasiSE
)
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.