Nothing
tests/testthat/test-make-outputs.R
In bage: Bayesian Estimation and Forecasting of Age-Specific Rates
## 'combine_stored_draws_point_inner_outer' -----------------------------------------
test_that("'combine_stored_draws_point_inner_outer' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:4,
sex = c("F", "M"),
region = c("a", "b"),
time = 2001:2005)
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + region * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, region:time ~ Lin())
mod <- set_n_draw(mod, n_draw = 10)
vars_inner <- c("age", "sex")
use_term <- make_use_term(mod = mod, vars_inner = vars_inner)
mod_inner <- reduce_model_terms(mod = mod, use_term = use_term)
mod_inner <- fit_default(mod_inner, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
mod_outer <- reduce_model_terms(mod = mod, use_term = !use_term)
mod_outer <- fit_default(mod_outer, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
mod_comb <- combine_stored_draws_point_inner_outer(mod = mod,
mod_inner = mod_inner,
mod_outer = mod_outer,
use_term = use_term)
terms <- make_terms_effectfree(mod)
is_inner <- terms %in% c("(Intercept)", "age", "sex", "age:sex")
is_outer <- terms %in% c("region", "time", "region:time")
expect_identical(mod_comb$draws_effectfree[is_inner, ],
mod_inner$draws_effectfree)
expect_identical(mod_comb$draws_effectfree[is_outer, ],
mod_outer$draws_effectfree)
expect_identical(ncol(mod_comb$draws_hyper),
ncol(mod_outer$draws_hyper))
expect_identical(ncol(mod_comb$draws_hyperrand),
ncol(mod_outer$draws_hyperrand))
expect_identical(mod_comb$point_effectfree[is_inner],
mod_inner$point_effectfree)
expect_identical(mod_comb$point_effectfree[is_outer],
mod_outer$point_effectfree)
})
## 'con_by_fitted' ------------------------------------------------------------
test_that("'con_by_fitted' works", {
set.seed(0)
prior <- RW()
fitted <- rvec::rnorm_rvec(n = 100, n_draw = 10)
along <- "time"
dimnames_term <- list(time = 2001:2010,
age = 0:4,
sex = 1:2)
var_time <- "time"
var_age <- "age"
ans <- con_by_fitted(prior = prior,
fitted = fitted,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age)
expect_equal(sum(ans[c(3, 13, 23, 33, 43)]),
sum(ans[c(4, 14, 24, 34, 44)]))
expect_equal(ans[10], -ans[60])
})
## 'draw_vals_components_fitted' ----------------------------------------------
test_that("'draw_vals_components_fitted' works", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp())
mod <- set_prior(mod, age:sex ~ SVD(HMD))
mod <- fit(mod)
set.seed(0)
ans <- draw_vals_components_fitted(mod)
expect_identical(names(ans), c("term", "component", "level", ".fitted"))
i <- ans$term == "age:sex" & ans$component == "effect"
expect_true(mean(abs(as.numeric(sum(ans$.fitted[i])))) > 0)
})
## 'generate_prior_helper' ----------------------------------------------------
test_that("'generate_prior_helper' works with valid inputs, not along-by", {
x <- N()
ans_obtained <- generate_prior_helper(x = x, n_element = 5, n_draw = 2)
draw <- factor(rep(c("Draw 1", "Draw 2"), each = 5))
ans_expected <- list(ans = tibble::tibble(draw = draw,
element = rep(1:5, times = 2)),
matrix_along_by = matrix(0:4, nc = 1, dimnames = list(seq_len(5), NULL)),
levels_effect = 1:5)
expect_identical(ans_obtained, ans_expected)
})
test_that("'generate_prior_helper' works with valid inputs, along-by, n_by = 1", {
x <- AR()
ans_obtained <- generate_prior_helper(x = x, n_along = 5, n_by = 1, n_draw = 2)
draw <- factor(rep(c("Draw 1", "Draw 2"), each = 5))
by <- factor(rep("By 1", 10))
levels_effect <- paste(1, 1:5, sep = ".")
ans_expected <- list(ans = tibble::tibble(draw = draw,
by = by,
along = rep(1:5, times = 2)),
matrix_along_by = matrix(0:4, nc = 1, dimnames = list(seq_len(5), 1L)),
levels_effect = levels_effect)
expect_identical(ans_obtained, ans_expected)
})
test_that("'generate_prior_helper' works with valid inputs, along-by, n_by = 2", {
x <- AR()
ans_obtained <- generate_prior_helper(x = x, n_along = 5, n_by = 2, n_draw = 2)
draw <- factor(rep(c("Draw 1", "Draw 2"), each = 10))
by <- factor(rep(rep(paste("By", 1:2), each = 5), times = 2))
ans_expected <- list(ans = tibble::tibble(draw = draw,
by = by,
along = rep(1:5, times = 4)),
matrix_along_by = matrix(0:9, nc = 2, dimnames = list(1:5, 1:2)),
levels_effect = paste(rep(1:2, each = 5), 1:5, sep = "."))
expect_identical(ans_obtained, ans_expected)
})
## 'generate_prior_svd_helper' ------------------------------------------------
test_that("'generate_prior_svd_helper' works with valid inputs - n_element", {
x <- SVD(LFP)
set.seed(0)
n_element <- 20
n_draw <- 25
ans_obtained <- generate_prior_svd_helper(x, n_element = n_element, n_draw = n_draw)
expect_true(is.list(ans_obtained))
})
test_that("'generate_prior_svd_helper' works with valid inputs - n_by = 1", {
x <- SVD_RW(LFP)
set.seed(0)
n_along <- 20
n_by <- 2
n_draw <- 25
ans_obtained <- generate_prior_svd_helper(x, n_along = n_along, n_by = n_by, n_draw = n_draw)
expect_true(is.list(ans_obtained))
})
## 'generate_ssvd_helper' -----------------------------------------------------------------
test_that("'generate_ssvd_helper' works with valid inputs - indep = NULL, n_element = 1", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_element = 1,
n_draw = 3,
n_comp = 2,
indep = NULL,
age_labels = NULL)
expect_identical(nrow(ans$matrix), length(unique(ans$ans$age)))
expect_identical(ncol(ans$matrix), 2L)
})
test_that("'generate_ssvd_helper' works with valid inputs - indep = TRUE, n_by = 1, n_along = 3", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_along = 3,
n_by = 1,
n_draw = 3,
n_comp = 2,
indep = TRUE,
age_labels = NULL)
expect_identical(nrow(ans$matrix), 3L * nrow(unique(ans$ans[c("age", "sexgender")])))
expect_identical(ncol(ans$matrix), 12L)
})
test_that("'generate_ssvd_helper' works with valid inputs - indep = FALSE, n_by = 1", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_along = 2,
n_by = 1,
n_draw = 3,
n_comp = 2,
indep = FALSE,
age_labels = NULL)
expect_identical(nrow(ans$matrix), 2L * nrow(unique(ans$ans[c("age", "sexgender")])))
expect_identical(ncol(ans$matrix), 4L)
})
test_that("'generate_ssvd_helper' works with valid inputs - indep = TRUE, n_by = 2", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_along = 3,
n_by = 2,
n_draw = 3,
n_comp = 2,
indep = TRUE,
age_labels = NULL)
expect_identical(nrow(ans$matrix), 3L * nrow(unique(ans$ans[c("by", "age", "sexgender")])))
expect_identical(ncol(ans$matrix), 24L)
expect_identical(nrow(ans$matrix_along_by), 3L)
})
## 'get_term_from_est' --------------------------------------------------------
test_that("'get_term_from_est' works", {
est <- list(effectfree = c(a = 1, a = 2, b = 3, c = 4, c = 5),
hyper = c(a = 1),
hyperrandfree = double(),
log_disp = c(disp = 3))
ans_obtained <- get_term_from_est(est = est, index_term = c(1L, 4L, 7L))
ans_expected <- c("a", "c", "disp")
expect_identical(ans_obtained, ans_expected)
})
## 'get_disp' -----------------------------------------------------------------
test_that("'get_disp' works - unfitted", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 5)
set.seed(1)
ans_obtained <- get_disp(mod)
set.seed(1)
ans_expected <- draw_vals_disp(mod, n_sim = 5)
expect_identical(ans_obtained, ans_expected)
})
test_that("'get_disp' works - fitted", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 5)
mod <- fit(mod)
set.seed(1)
ans_obtained <- get_disp(mod)
set.seed(1)
ans_expected <- rvec::rvec(matrix(mod$draws_disp, nr = 1))
expect_identical(ans_obtained, ans_expected)
})
## 'insert_after' -------------------------------------------------------
test_that("'insert_after' works with data frames", {
df <- data.frame(x = 1:3, y = 3:1)
x <- 11:13
nm_x = "new"
ans_obtained <- insert_after(df = df,
nm_after = "x",
x = x,
nm_x = nm_x)
ans_expected <- data.frame(x = 1:3, new = 11:13, y = 3:1)
expect_identical(ans_obtained, ans_expected)
ans_obtained <- insert_after(df = df,
nm_after = "y",
x = x,
nm_x = nm_x)
ans_expected <- data.frame(x = 1:3, y = 3:1, new = 11:13)
expect_identical(ans_obtained, ans_expected)
})
test_that("'insert_after' works with tibbles", {
df <- tibble::tibble(x = 1:3, y = 3:1)
x <- 11:13
nm_x = "new"
ans_obtained <- insert_after(df = df,
nm_after = "x",
x = x,
nm_x = nm_x)
ans_expected <- tibble(x = 1:3, new = 11:13, y = 3:1)
expect_identical(ans_obtained, ans_expected)
ans_obtained <- insert_after(df = df,
nm_after = "y",
x = x,
nm_x = nm_x)
ans_expected <- tibble(x = 1:3, y = 3:1, new = 11:13)
expect_identical(ans_obtained, ans_expected)
})
## 'is_same_class' ------------------------------------------------------------
test_that("'is_same_class' returns TRUE when classes same", {
expect_true(is_same_class(AR1(), AR1()))
expect_true(is_same_class(1L, 2L))
})
test_that("'is_same_class' returns FALSE when classes different", {
expect_false(is_same_class(AR1(), N()))
expect_false(is_same_class(1L, FALSE))
})
## 'make_combined_matrix_effect_outcome' -----------------------------------------
test_that("'make_combined_matrix_effect_outcome' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_combined_matrix_effect_outcome(mod)
expect_identical(nrow(ans_obtained), nrow(data))
expect_identical(ncol(ans_obtained), length(make_terms_effects(mod$dimnames_terms)))
expect_false(all(ans_obtained[,1] == 0))
})
## 'make_combined_matrix_effectfree_effect' -----------------------------------------
test_that("'make_combined_matrix_effectfree_effect' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_combined_matrix_effectfree_effect(mod)
expect_identical(nrow(ans_obtained), length(make_terms_effects(mod$dimnames_terms)))
expect_identical(ncol(ans_obtained), length(make_terms_effectfree(mod)))
})
## 'make_combined_offset_effectfree_effect' -----------------------------------------
test_that("'make_combined_offset_effectfree_effect' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_combined_offset_effectfree_effect(mod)
expect_identical(length(ans_obtained), length(make_terms_effects(mod$dimnames_terms)))
expect_true(all(ans_obtained == 0))
})
## 'make_comp_components' -----------------------------------------------------
test_that("'make_comp_components' works - no hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
draws <- make_draws_components(mod)
ans <- make_comp_components(mod)
expect_identical(length(draws), length(ans))
expect_setequal(ans, c("effect", "hyper", "disp"))
})
test_that("'make_comp_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ Lin())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
draws <- make_draws_components(mod)
ans <- make_comp_components(mod)
expect_identical(length(draws), length(ans))
expect_setequal(ans, c("effect", "hyper", "trend", "error", "disp"))
})
test_that("'make_comp_components' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD))
mod <- set_prior(mod, age:sex ~ Sp())
set.seed(0)
mod <- fit(mod)
term <- make_term_components(mod)
ans <- make_comp_components(mod)
expect_identical(length(term), length(ans))
expect_setequal(ans, c("effect", "hyper", "svd", "spline", "disp"))
})
## 'make_comp_hyperrand' ------------------------------------------------------
test_that("'make_comp_hyperrand' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex:time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
ans_obtained <- make_comp_hyperrand(mod)
ans_expected <- c("hyper", "hyper", rep(c("trend", "error"), each = 12))
expect_identical(ans_obtained, ans_expected)
})
## 'make_copies_repdata' ------------------------------------------------------
test_that("'make_copies_repdata' works with valid inputs", {
data <- tibble::tibble(age = 0:4, deaths = 1:5, population = 11:15)
n <- 2
ans_obtained <- make_copies_repdata(data = data, n = n)
ans_expected <- tibble::tibble(.replicate = factor(rep(c("Original", "Replicate 1", "Replicate 2"),
each = 5),
levels = c("Original",
"Replicate 1", "Replicate 2")),
age = rep(0:4, 3),
deaths = rep(1:5, 3),
population = rep(11:15, 3))
expect_identical(ans_obtained, ans_expected)
})
## 'make_draws_components' -----------------------------------------------
test_that("'make_draws_components' works - no svd, spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
expect_identical(length(ans_obtained),
length(make_effectfree(mod)) +
length(make_hyper(mod)) +
length(make_hyperrand(mod)) + 1L)
})
test_that("'make_draws_components' works - has spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 6))
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
expect_identical(length(ans_obtained),
length(make_effectfree(mod)) + length(unique(data$age)) +
length(make_hyper(mod)) +
length(make_hyperrand(mod)) + 1L)
})
test_that("'make_draws_components' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD, n_comp = 3))
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
expect_identical(length(ans_obtained),
length(make_effectfree(mod)) + length(unique(data$age)) +
length(make_hyper(mod)) +
length(make_hyperrand(mod)) + 1L)
})
test_that("'make_draws_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age:sex ~ Lin())
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
})
## 'make_draws_disp' ----------------------------------------------------
test_that("'make_draws_disp' works", {
set.seed(0)
draws_post <- matrix(rnorm(13 * 5), nrow = 13)
ans_obtained <- make_draws_disp(draws_post)
ans_expected <- exp(draws_post[13,])
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_draws_effectfree' ----------------------------------------------------
test_that("'make_draws_effectfree' works", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = rnorm(2),
hyperrand = numeric(),
log_disp = 0.5)
draws_post <- matrix(rnorm(13 * 5), nrow = 13)
ans_obtained <- make_draws_effectfree(est = est, draws_post = draws_post)
ans_expected <- draws_post[1:10,]
expect_identical(ans_obtained, ans_expected)
})
## 'make_draws_hyper' ----------------------------------------------------
test_that("'make_draws_hyper' works", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = runif(4),
hyperrand = numeric(),
log_disp = 0.5)
transforms_hyper <- list(identity, identity, exp, exp)
draws_post <- matrix(rnorm(15 * 5), nrow = 15)
ans_obtained <- make_draws_hyper(est = est,
transforms_hyper = transforms_hyper,
draws_post = draws_post)
ans_expected <- draws_post[11:14, ]
ans_expected[3:4,] <- exp(ans_expected[3:4,])
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_draws_hyperrandfree' -------------------------------------------------
test_that("'make_draws_hyperrandfree' works - has hyperrandfree", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = runif(4),
hyperrandfree = rnorm(5),
log_disp = 0.5)
draws_post <- matrix(rnorm(20 * 5), nrow = 20)
ans_obtained <- make_draws_hyperrandfree(est = est,
draws_post = draws_post)
ans_expected <- draws_post[15:19, ]
expect_identical(unname(ans_obtained), ans_expected)
})
test_that("'make_draws_hyperrandfree' works - no hyperrandfree", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = runif(4),
hyperrandfree = numeric(),
log_disp = 0.5)
draws_post <- matrix(rnorm(15 * 5), nrow = 15)
ans_obtained <- make_draws_hyperrandfree(est = est,
draws_post = draws_post)
ans_expected <- matrix(0, nrow = 0, ncol = 5)
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_draws_post' ------------------------------------------------------
test_that("'make_draws_post' works with valid inputs - has R_prec", {
set.seed(0)
est <- list(effectfree = c("(Intercept)" = -3,
sex = c(-1, 1),
age = rnorm(10)),
log_disp = 0.2)
prec <- diag(12)
prec <- Matrix::Matrix(prec)
map <- list(effectfree = factor(c(1, NA, NA, 2:11)),
disp = factor(1))
ans <- make_draws_post(est = est,
prec = prec,
map = map,
n_draw = 500000L)
expect_equal(rowMeans(ans), unlist(est, use.names = FALSE), tolerance = 0.02)
expect_equal(solve(cov(t(ans[c(1, 4:14),]))), as.matrix(prec), tolerance = 0.02)
})
## 'fit_default' --------------------------------------------------------------
## lots more tests for 'fit' method
test_that("'fit_default' works with pois, optimzier is 'multi'", {
set.seed(10)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(age ~ AR()) |>
set_prior(age:sex ~ RW2(sd = 0)) |>
set_prior(age:sex:time ~ AR())
ans_obtained <- fit_default(mod, optimizer = "multi", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
expect_s3_class(ans_obtained, "bage_mod")
})
test_that("'fit_default' works with pois, optimzier is 'nlminb'", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
expect_s3_class(ans_obtained, "bage_mod")
})
test_that("'fit_default' works with pois - start_oldpar", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
ans_obtained <- fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = TRUE)
expect_s3_class(ans_obtained, "bage_mod")
})
test_that("'fit_default' gives error with 'start_oldpar' if model fitted", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
expect_error(fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = TRUE),
"`start_oldpar` is TRUE but model has not been fitted.")
})
## 'fit_inner_outer' ----------------------------------------------------------
test_that("'fit_inner_outer' works with with pois", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$popn <- rpois(n = nrow(data), lambda = 1000)
data$deaths <- rnbinom(n = nrow(data), mu = 0.1 * data$popn, size = 100)
formula <- deaths ~ age * sex + region * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
set.seed(0)
ans_inner_outer <- fit_inner_outer(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
vars_inner = c("age", "sex"))
set.seed(0)
ans_default <- fit_default(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
aggregate = TRUE)
aug_inner_outer <- ans_inner_outer |>
augment()
fit_inner_outer <- rvec::draws_median(aug_inner_outer$.fitted)
aug_default <- ans_default |>
augment()
fit_default <- rvec::draws_median(aug_default$.fitted)
expect_true(cor(fit_inner_outer, fit_default) > 0.98)
})
test_that("'fit_inner_outer' works with with binom", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$popn <- rpois(n = nrow(data), lambda = 1000)
data$deaths <- rbinom(n = nrow(data), prob = 0.2, size = data$popn)
formula <- deaths ~ age * sex + region * time
mod <- mod_binom(formula = formula,
data = data,
size = popn)
set.seed(0)
ans_inner_outer <- fit_inner_outer(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
vars_inner = NULL)
set.seed(0)
ans_default <- fit_default(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
aggregate = TRUE)
aug_inner_outer <- ans_inner_outer |>
augment()
fit_inner_outer <- rvec::draws_median(aug_inner_outer$.fitted)
aug_default <- ans_default |>
augment()
fit_default <- rvec::draws_median(aug_default$.fitted)
expect_true(cor(fit_inner_outer, fit_default) > 0.98)
})
test_that("'fit_inner_outer' works with with norm", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$wt <- rpois(n = nrow(data), lambda = 10)
data$income <- rnorm(n = nrow(data), mean = data$age + data$time/100, sd = 5 / sqrt(data$wt))
formula <- income ~ age * sex + region * time
mod <- mod_norm(formula = formula,
data = data,
weights = wt)
set.seed(0)
ans_inner_outer <- fit_inner_outer(mod,
optimizer = "BFGS",
quiet = TRUE,
start_oldpar = FALSE,
vars_inner = c("age", "sex"))
set.seed(0)
ans_default <- fit_default(mod,
optimizer = "BFGS",
quiet = TRUE,
start_oldpar = FALSE,
aggregate = TRUE)
aug_inner_outer <- ans_inner_outer |>
augment()
fit_inner_outer <- rvec::draws_median(aug_inner_outer$.fitted)
aug_default <- ans_default |>
augment()
fit_default <- rvec::draws_median(aug_default$.fitted)
expect_true(cor(fit_inner_outer, fit_default) > 0.98)
})
test_that("'fit_inner_outer' throws error when 'start_oldpar' is TRUE", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$wt <- rpois(n = nrow(data), lambda = 10)
data$income <- rnorm(n = nrow(data), mean = data$age + data$time/100, sd = 5 / sqrt(data$wt))
formula <- income ~ age * sex + region * time
mod <- mod_norm(formula = formula,
data = data,
weights = wt)
set.seed(0)
expect_error(fit_inner_outer(mod,
optimizer = "BFGS",
quiet = TRUE,
start_oldpar = TRUE,
vars_inner = c("age", "sex")),
"`start_oldpar` must be FALSE when using \"inner-outer\" method.")
})
## 'make_along_mod' -----------------------------------------------------------
test_that("'make_along_mod' works", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:time ~ Sp(n_comp = 5, along = "age"))
ans_obtained <- make_along_mod(mod)
ans_expected <- c("(Intercept)" = NA,
sex = NA,
age = "age",
time = "time",
"sex:age" = "age",
"age:time" = "age")
expect_identical(ans_obtained, ans_expected)
})
## 'make_hyperrand' -----------------------------------------------------------
test_that("'make_hyperrand' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ Lin())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
ans <- make_hyperrand(mod)
expect_identical(length(ans), 26L)
})
test_that("'make_hyperrand' works - no hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ AR())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
ans <- make_hyperrand(mod)
expect_identical(length(ans), 0L)
})
## 'make_hyperrand_lin' -------------------------------------------------------
test_that("'make_hyperrand_lin' works with main effect", {
prior <- Lin()
hyperrandfree <- rvec::rnorm_rvec(n = 1, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_lin(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
trend <- hyperrandfree * (1:10 - mean(1:10))
error <- effectfree - trend
slope <- trend[2] - trend[1]
ans_expected <- c(slope, trend, error)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_lin' works with interaction, con is 'none'", {
prior <- Lin_AR()
hyperrandfree <- rvec::rnorm_rvec(n = 2, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_lin(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
trend <- c(hyperrandfree[1] * (1:10 - mean(1:10)),
hyperrandfree[2] * (1:10 - mean(1:10)))
error <- effectfree - trend
slope <- c(trend[2] - trend[1], trend[12] - trend[11])
ans_expected <- c(slope, trend, error)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_lin' works with interaction, con is 'by'", {
prior <- Lin_AR(con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 1, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_lin(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
trend <- c(-sqrt(0.5) * hyperrandfree * (1:10 - mean(1:10)),
sqrt(0.5) * hyperrandfree * (1:10 - mean(1:10)))
error <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - trend
slope <- c(trend[2] - trend[1], trend[12] - trend[11])
ans_expected <- c(slope, trend, error)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_randomseasfix' -----------------------------------------------
test_that("'make_hyperrand_randomseasfix' works with main effect", {
prior <- RW_Seas(n = 4)
hyperrandfree <- rvec::rnorm_rvec(n = 3, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(hyperrandfree[1],
hyperrandfree[2],
hyperrandfree[3],
-sum(hyperrandfree))[c(1:4, 1:4, 1:2)]
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasfix' works with interaction, random_sum is FALSE", {
set.seed(0)
prior <- RW2_Seas(n = 4)
hyperrandfree <- rvec::rnorm_rvec(n = 6, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(hyperrandfree[1],
hyperrandfree[2],
hyperrandfree[3],
-sum(hyperrandfree[1:3]))[c(1:4, 1:4, 1:2)]
season2 <- c(hyperrandfree[4],
hyperrandfree[5],
hyperrandfree[6],
-sum(hyperrandfree[4:6]))[c(1:4, 1:4, 1:2)]
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasfix' works with interaction, con is 'by'", {
prior <- RW_Seas(n = 4, con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 3, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(hyperrandfree[1],
hyperrandfree[2],
hyperrandfree[3],
-sum(hyperrandfree[1:3]))[c(1:4, 1:4, 1:2)]
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_randomseasvary' ----------------------------------------------
test_that("'make_hyperrand_randomseasvary' works with main effect", {
prior <- RW_Seas(n = 4, s_seas = 1)
hyperrandfree <- rvec::rnorm_rvec(n = 8, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(hyperrandfree[1:3],
-sum(hyperrandfree[1:3]),
hyperrandfree[4:6],
-sum(hyperrandfree[4:6]),
hyperrandfree[7:8])
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasvary' works with interaction, con is 'none'", {
set.seed(0)
prior <- RW2_Seas(n = 4, s_seas = 1)
hyperrandfree <- rvec::rnorm_rvec(n = 16, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(hyperrandfree[1:3],
-sum(hyperrandfree[1:3]),
hyperrandfree[4:6],
-sum(hyperrandfree[4:6]),
hyperrandfree[7:8])
season2 <- c(hyperrandfree[9:11],
-sum(hyperrandfree[9:11]),
hyperrandfree[12:14],
-sum(hyperrandfree[12:14]),
hyperrandfree[15:16])
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasvary' works with interaction, con is 'by'", {
prior <- RW_Seas(n_seas = 4, s_seas = 1, con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 8, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(hyperrandfree[1:3],
-sum(hyperrandfree[1:3]),
hyperrandfree[4:6],
-sum(hyperrandfree[4:6]),
hyperrandfree[7:8])
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_zeroseasfix' -----------------------------------------------
test_that("'make_hyperrand_zeroseasfix' works with main effect", {
prior <- RW_Seas(n = 4, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 2, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree))[c(1:4, 1:4, 1:2)]
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasfix' works with interaction, con is 'none'", {
set.seed(0)
prior <- RW2_Seas(n = 4, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 4, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]))[c(1:4, 1:4, 1:2)]
season2 <- c(effectfree[11],
effectfree[11] + hyperrandfree[3],
effectfree[11] + hyperrandfree[4],
-3 * effectfree[11] - sum(hyperrandfree[3:4]))[c(1:4, 1:4, 1:2)]
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasfix' works with interaction, con is 'by'", {
prior <- RW_Seas(n = 4, sd = 0, con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 2, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]))[c(1:4, 1:4, 1:2)]
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_zeroseasvary' ----------------------------------------------
test_that("'make_hyperrand_zeroseasvary' works with main effect", {
prior <- RW_Seas(n = 4, s_seas = 1, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 7, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]),
effectfree[1] + hyperrandfree[3:5],
-3 * effectfree[1] - sum(hyperrandfree[3:5]),
effectfree[1] + hyperrandfree[6:7])
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasvary' works with interaction, con is 'none'", {
set.seed(0)
prior <- RW2_Seas(n = 4, s_seas = 1, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 14, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]),
effectfree[1] + hyperrandfree[3:5],
-3 * effectfree[1] - sum(hyperrandfree[3:5]),
effectfree[1] + hyperrandfree[6:7])
season2 <- c(effectfree[11],
effectfree[11] + hyperrandfree[8],
effectfree[11] + hyperrandfree[9],
-3 * effectfree[11] - sum(hyperrandfree[8:9]),
effectfree[11] + hyperrandfree[10:12],
-3 * effectfree[11] - sum(hyperrandfree[10:12]),
effectfree[11] + hyperrandfree[13:14])
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasvary' works with interaction, con is 'by'", {
prior <- RW_Seas(n_seas = 4, s_seas = 1, con = "by", sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 7, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]),
effectfree[1] + hyperrandfree[3:5],
-3 * effectfree[1] - sum(hyperrandfree[3:5]),
effectfree[1] + hyperrandfree[6:7])
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_levels_spline' -------------------------------------------------------
test_that("'make_levels_spline' works - unlist is FALSE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:sex ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:time ~ Sp(n_comp = 5, along = "age"))
set.seed(0)
ans_obtained <- make_levels_spline(mod, unlist = FALSE)
ans_expected <- list("(Intercept)" = NULL,
sex = NULL,
age = paste0("comp", 1:5),
time = NULL,
"sex:age" = paste(c("F", "M"), rep(paste0("comp", 1:5), each = 2), sep = "."),
"age:time" = paste0("comp", 1:5, ".", rep(2000:2005, each = 5)))
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_spline' works - unlist is TRUE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:sex ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:time ~ Sp(n_comp = 5, along = "age"))
set.seed(0)
ans_obtained <- make_levels_spline(mod, unlist = TRUE)
ans_expected <- c(paste0("comp", 1:5),
paste(c("F", "M"), rep(paste0("comp", 1:5), each = 2), sep = "."),
paste(paste0("comp", 1:5), rep(2000:2005, each = 5), sep = "."))
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_spline_term' --------------------------------------------------
test_that("'make_levels_spline_term' works - con is 'none'", {
prior <- Sp(n = 5)
dimnames_term <- list(reg = 1:2,
age = 1:20)
ans_obtained <- make_levels_spline_term(prior = prior,
dimnames_term = dimnames_term,
var_time = "time",
var_age = "age")
ans_expected <- paste(1:2, rep(paste0("comp", 1:5), each = 2), sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_spline_term' works - con is 'by'", {
prior <- Sp(n = 5, con = "by")
dimnames_term <- list(reg = 1:2,
age = 1:20)
ans_obtained <- make_levels_spline_term(prior = prior,
dimnames_term = dimnames_term,
var_time = "time",
var_age = "age")
ans_expected <- paste("reg1", paste0("comp", 1:5), sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_spline_term' works - con is 'by', more dimensions", {
prior <- Sp(n = 5, con = "by")
dimnames_term <- list(reg = 1:4,
age = 1:20,
sex = c("f", "m"))
ans_obtained <- make_levels_spline_term(prior = prior,
dimnames_term = dimnames_term,
var_time = "time",
var_age = "age")
ans_expected <- paste(paste0("reg", 1:3),
rep(paste0("comp", 1:5), each = 3),
rep("sex1", times = 15),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_svd' ----------------------------------------------------------
test_that("'make_levels_svd' works - unlist is FALSE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD))
mod <- set_prior(mod, age:sex ~ SVD(HMD))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD))
set.seed(0)
ans_obtained <- make_levels_svd(mod, unlist = FALSE)
ans_expected <- list("(Intercept)" = NULL,
sex = NULL,
age = paste0("comp", 1:3),
time = NULL,
"sex:age" = paste0(rep(c("F", "M"), each = 3), ".comp", 1:3),
"age:time" = paste0("comp", 1:3, ".", rep(2000:2005, each = 3)))
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd' works - unlist is TRUE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:sex ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD, n_comp = 5))
set.seed(0)
ans_obtained <- make_levels_svd(mod, unlist = TRUE)
ans_expected <- c(paste0("comp", 1:5),
paste0(rep(c("F", "M"), each = 5), ".comp", 1:5),
paste0("comp", 1:5, ".", rep(2000:2005, each = 5)))
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_svd_term' ----------------------------------------------------------
test_that("'make_levels_svd_term' works - total, main effect", {
prior <- SVD(HMD)
dimnames_term <- list(age = c(0:59, "60+"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste0("comp", 1:3)
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - joint, age:sex", {
prior <- SVD(HMD, n_comp = 5)
dimnames_term <- list(age = c(0:59, "60+"), sex = c("M", "F"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(rep(c("M", "F"), each = 5), paste0("comp", 1:5), sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - indep, age:sex:reg", {
prior <- SVD(HMD, indep = FALSE, n_comp = 5)
dimnames_term <- list(reg = 1:2, age = c(0:59, "60+"), sex = c("M", "F"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(paste0("comp", 1:5),
rep(1:2, each = 5),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - indep, age:sex:time, con is 'none'", {
prior <- SVD_RW(HMD, indep = FALSE, n_comp = 5)
dimnames_term <- list(time = 2001:2010, age = c(0:59, "60+"), sex = c("M", "F"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(paste0("comp", 1:5),
rep(2001:2010, each = 5),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - indep, age:time:reg, con is 'by'", {
prior <- SVD_RW(HMD, indep = FALSE, n_comp = 5, con = "by")
dimnames_term <- list(time = 2001:2010, age = c(0:59, "60+"), reg = c("a", "b", "c"))
var_age <- "age"
var_sexgender <- "sex"
var_time <- "time"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(paste0("comp", 1:5),
rep(2001:2010, each = 5),
rep(c("reg1", "reg2"), each = 50),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
## 'make_lin_trend' -----------------------------------------------------------
test_that("'make_lin_trend' works with valid inputs - n_by = 1", {
slope <- rvec::rvec(matrix(as.numeric(1:5), nr = 1))
matrix_along_by <- matrix(0:9, nr = 10)
ans_obtained <- make_lin_trend(slope = slope,
matrix_along_by = matrix_along_by)
ones <- rep(1, times = 10)
s <- 1:10
ans_expected <- rvec::rvec(outer(ones, -0.5 * 11 * (1:5)) + outer(s, 1:5))
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_lin_trend' works with valid inputs - n_by = 2, transposed", {
slope <- rvec::rvec(matrix(rnorm(10), nr = 2))
intercept <- -0.5 * 6 * slope
matrix_along_by <- t(matrix(0:9, nr = 2))
ans_obtained <- make_lin_trend(slope = slope,
matrix_along_by = matrix_along_by)
s <- 1:5
ans_expected <- c(intercept[1] + slope[1] * s, intercept[2] + slope[2] * s)
ans_expected <- ans_expected[c(1, 6, 2, 7, 3, 8, 4, 9, 5, 10)]
expect_identical(ans_obtained, ans_expected)
})
## 'make_stored_draws' --------------------------------------------------------
test_that("'make_stored_draws' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex ~ Known(c(0.1, -0.1)))
mod <- set_n_draw(mod, n = 10)
est <- list(effectfree = c(rnorm(11), 0.1, -0.1),
hyper = rnorm(1),
hyperrandfree = numeric(),
disp = runif(1))
prec <- crossprod(matrix(rnorm(169), nr = 13))
map <- make_map(mod)
ans <- make_stored_draws(mod = mod,
est = est,
prec = prec,
map = map)
expect_identical(ncol(ans$draws_effectfree), 10L)
expect_identical(ncol(ans$draws_hyper), 10L)
expect_identical(length(ans$draws_disp), 10L)
ans2 <- make_stored_draws(mod = mod,
est = est,
prec = prec,
map = map)
expect_identical(ans, ans2)
})
## 'make_stored_point' --------------------------------------------------------
test_that("'make_stored_draws' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex ~ Known(c(0.1, -0.1)))
est <- list(effectfree = c(rnorm(11), 0.1, -0.1),
hyper = rnorm(1),
hyperrandfree = numeric(),
log_disp = runif(1))
ans <- make_stored_point(mod = mod,
est = est)
expect_identical(ans$point_effectfree, est$effectfree)
expect_identical(ans$point_hyper, exp(est$hyper))
expect_identical(ans$point_hyperrandfree, double())
expect_identical(ans$point_disp, exp(est$log_disp))
})
## 'make_par_disp' ------------------------------------------------------------
test_that("'make_par_disp' works with bage_mod_pois", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + time + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 10)
mod <- fit(mod)
components <- components(mod)
meanpar <- exp(make_linpred_comp(components = components,
data = mod$data,
dimnames_terms = mod$dimnames_terms))
disp <- components$.fitted[components$component == "disp"]
set.seed(1)
ans_obtained <- make_par_disp(mod,
meanpar = meanpar,
disp = disp)
set.seed(1)
ans_expected <- rvec::rgamma_rvec(n = length(meanpar),
data$deaths + 1/disp,
data$popn + 1/(disp*meanpar))
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_par_disp' works with bage_mod_binom", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rbinom(n = nrow(data), size = data$popn, prob = 0.3)
formula <- deaths ~ age + time + sex
mod <- mod_binom(formula = formula,
data = data,
size = popn)
mod <- fit(mod)
mod <- set_n_draw(mod, n = 10)
components <- components(mod)
invlogit <- function(x) 1 / (1 + exp(-x))
meanpar <- invlogit(make_linpred_comp(components = components,
data = mod$data,
dimnames_terms = mod$dimnames_terms))
disp <- components$.fitted[components$component == "disp"]
set.seed(1)
ans_obtained <- make_par_disp(mod,
meanpar = meanpar,
disp = disp)
set.seed(1)
ans_expected <- rvec::rbeta_rvec(n = length(meanpar),
data$deaths + meanpar/disp,
data$popn - data$deaths + (1 - meanpar)/disp)
expect_equal(ans_obtained, ans_expected)
})
## 'make_is_fixed' ------------------------------------------------------------
test_that("'make_is_fixed' works when nothing fixed", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- fit(mod)
est <- list(effectfree = make_effectfree(mod),
hyper = make_hyper(mod),
hyperrand = make_hyperrand(mod),
log_disp = 0)
map <- make_map(mod)
expect_true(is.null(map))
ans_obtained <- make_is_fixed(est = est, map = map)
ans_expected <- rep(FALSE, times = length(unlist(est)))
})
test_that("'make_is_fixed' works when Known prior", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex ~ Known(c(0.1, -0.1)))
mod <- fit(mod)
est <- list(effectfree = make_effectfree(mod),
hyper = make_hyper(mod),
hyperrand = make_hyperrand(mod),
log_disp = 0)
map <- make_map(mod)
ans_obtained <- make_is_fixed(est = est, map = map)
ans_expected <- rep(c(FALSE, TRUE, FALSE),
times = c(11,
2,
20 + 6 + length(est$hyper) +
+ length(est$log_disp)))
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_level_components' ----------------------------------------------------
test_that("'make_level_components' works - no hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_disp(mod, mean = 0)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_level_components(mod)
expect_identical(length(ans), length(comp))
})
test_that("'make_level_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_disp(mod, mean = 0)
mod <- set_prior(mod, sex:time ~ Lin())
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_level_components(mod)
expect_identical(length(ans), length(comp))
})
## 'make_levels_hyper' --------------------------------------------------------
test_that("'make_levels_hyper' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_levels_hyper(mod)
ans_expected <- c(age = "sd",
time = "sd",
"age:sex" = "sd")
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_hyperrand' ----------------------------------------------------
test_that("'make_levels_hyperrand' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex:time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
ans_obtained <- make_levels_hyperrand(mod, unlist = TRUE)
ans_expected <- c("slope.F", "slope.M",
rep(paste(c("F", "M"), rep(2000:2005, each = 2), sep = "."), 2))
expect_identical(ans_obtained, ans_expected)
ans_obtained <- make_levels_hyperrand(mod, unlist = FALSE)
ans_expected <- list(character(),
character(),
c("slope.F", "slope.M",
rep(paste(c("F", "M"), rep(2000:2005, each = 2), sep = "."), 2)))
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_replicate' ----------------------------------------------------
test_that("'make_levels_replicate' works", {
ans_obtained <- make_levels_replicate(n = 2, n_row_data = 3)
ans_expected <- factor(rep(c("Original", "Replicate 1", "Replicate 2"),
each = 3),
levels = c("Original", "Replicate 1", "Replicate 2"))
expect_identical(ans_obtained, ans_expected)
})
## 'make_linpred_comp' --------------------------------------------------------
test_that("'make_linpred_comp' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 10L)
comp <- components(mod, quiet = TRUE)
ans <- make_linpred_comp(components = comp,
data = mod$data,
dimnames_terms = mod$dimnames_terms)
expect_identical(length(ans), length(mod$outcome))
})
## 'make_linpred_raw' ---------------------------------------------------------
test_that("'make_linpred_raw' works with valid inputs - point is FALSE", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 10L)
mod <- fit(mod)
ans_obtained <- make_linpred_raw(mod, point = FALSE)
comp <- components(mod, quiet = TRUE)
ans_expected <- make_linpred_comp(components = comp,
data = mod$data,
dimnames_terms = mod$dimnames_terms)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_linpred_raw' works with valid inputs - point is TRUE", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 10L)
mod <- fit(mod)
ans_obtained <- make_linpred_raw(mod, point = TRUE)
m1 <- make_combined_matrix_effect_outcome(mod)
m2 <- make_combined_matrix_effectfree_effect(mod)
ans_expected <- as.double(m1 %*% m2 %*% mod$point_effectfree)
expect_equal(ans_obtained, ans_expected)
})
## 'make_point_est_effects' ---------------------------------------------------
test_that("'make_point_est_effects' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(age ~ RW(sd = 0)) |>
set_prior(age:sex ~ RW(sd = 0))
mod <- fit(mod)
ans_obtained <- make_point_est_effects(mod)
int <- unname(mod$point_effectfree[1])
age <- c(0, unname(mod$point_effectfree[2:10]))
sex <- unname(mod$point_effectfree[11:12])
agesex <- c(0, unname(mod$point_effectfree[13:21]),
0, unname(mod$point_effectfree[22:30]))
ans_expected <- list("(Intercept)" = int,
age = age,
sex = sex,
"age:sex" = agesex)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_point_est_effects' throws correct error when not fitted", {
set.seed(0)
data <- expand.grid(age = 0:9,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
expect_error(make_point_est_effects(mod),
"Internal error: Model not fitted.")
})
## 'make_scaled_eigen' --------------------------------------------------------
## See also tests for rvnorm_eigen
test_that("'make_scaled_eigen' works with positive definite matrix", {
set.seed(0)
prec <- solve(crossprod(matrix(rnorm(25), 5)))
ans <- make_scaled_eigen(prec)
expect_identical(dim(ans), dim(prec))
})
test_that("'make_scaled_eigen' works with non-negative definite matrix", {
set.seed(0)
prec <- solve(crossprod(matrix(rnorm(25), 5)))
prec[5,] <- 0
prec[,5] <- 0
ans <- make_scaled_eigen(prec)
expect_identical(dim(ans), dim(prec))
})
## 'make_spline' --------------------------------------------------------------
test_that("'make_spline' works", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age ~ RW(sd = 0))
mod <- set_prior(mod, age:sex ~ Sp(n = 5))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_spline(mod = mod, effectfree = effectfree)
ans_expected <- effectfree[23:32,]
expect_equal(ans_obtained, ans_expected)
})
## 'make_svd' -----------------------------------------------------------------
test_that("'make_svd' works - SVD_RW, random", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:sex ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD, n_comp = 5))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_svd(mod = mod, effectfree = effectfree)
ans_expected <- effectfree[24:63,]
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_svd' works - SVD_RW, zero", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:sex ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD, n_comp = 5, sd = 0))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_svd(mod = mod, effectfree = effectfree)
ans_expected <- rbind(effectfree[24:33,],
matrix(0, nrow = 5, ncol = 5),
effectfree[34:58,])
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_svd' works - SVD_AR", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 20)
formula <- deaths ~ age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:time ~ SVD_AR1(HMD, n_comp = 1))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_svd(mod = mod, effectfree = effectfree)
ans_expected <- effectfree[22:27,]
expect_equal(ans_obtained, ans_expected)
})
## 'make_term_components' -----------------------------------------------------
test_that("'make_term_components' works - no disp", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_disp(mod, mean = 0)
mod <- set_prior(mod, age ~ Sp())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_term_components(mod)
expect_identical(length(ans), length(comp))
})
test_that("'make_term_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_term_components(mod)
expect_identical(length(ans), length(comp))
})
test_that("'make_term_components' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_term_components(mod)
expect_identical(length(ans), length(comp))
})
## 'make_term_spline' ------------------------------------------------------------
test_that("'make_term_spline' works - no spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
ans_obtained <- make_term_spline(mod)
ans_expected <- factor()
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_term_spline' works - has spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 5))
ans_obtained <- make_term_spline(mod)
ans_expected <- factor(rep("age", times = 5))
expect_identical(ans_obtained, ans_expected)
})
## 'make_term_svd' ------------------------------------------------------------
test_that("'make_term_svd' works - no svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 3)
ans_obtained <- make_term_svd(mod)
ans_expected <- factor()
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_term_svd' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD))
ans_obtained <- make_term_svd(mod)
ans_expected <- factor(c(rep("age", times = 3),
rep("age:time", times = 3 * 6)))
expect_identical(ans_obtained, ans_expected)
})
## 'make_terms_hyperrand' ----------------------------------------------------
test_that("'make_terms_hyperrand' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex:time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
ans_obtained <- make_terms_hyperrand(mod)
ans_expected <- factor(rep("sex:time", 26),
levels = c("(Intercept)", "age", "sex:time"))
expect_identical(ans_obtained, ans_expected)
})
## 'make_unconstr_dimnames_by' ------------------------------------------------
test_that("'make_unconstr_dimnames_by' works", {
dimnames_term <- list(age = 1:5, time = 1:6, reg = 1:3)
ans_obtained <- make_unconstr_dimnames_by(i_along = 2L,
dimnames_term = dimnames_term)
ans_expected <- list(age = paste0("age", 1:4), reg = paste0("reg", 1:2))
expect_identical(ans_obtained, ans_expected)
})
## 'paste_dot' ----------------------------------------------------------------
test_that("'paste_dot' works with valid inputs", {
expect_identical(paste_dot(1:3, 3:1), c("1.3", "2.2", "3.1"))
})
## 'rescale_lin_intercept' ----------------------------------------------------
test_that("'rescale_lin_intercept' works with valid inputs", {
set.seed(0)
slope <- rvec::rnorm_rvec(n = 4, n_draw = 10)
effect <- rvec::rnorm_rvec(n = 20, n_draw = 10)
matrix_along_by <- matrix(0:19, nrow = 5)
ans_obtained <- rescale_lin_intercept(slope = slope,
effect = effect,
matrix_along_by = matrix_along_by)
ans_expected <- slope
for (i in 1:4)
ans_expected[i] <- mean(effect[1:5 + (i - 1) * 5]) - mean(slope[i] * (1:5))
expect_equal(ans_obtained, ans_expected)
})
## 'make_transforms_hyper' ----------------------------------------------------
test_that("'make_transforms_hyper' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- fit(mod)
ans_obtained <- make_transforms_hyper(mod)
ans_expected <- list(exp, exp)
expect_identical(unname(ans_obtained), ans_expected,
ignore_function_env = TRUE)
})
## 'infer_trend_cyc_seas_err_forecast' --------------------------------------------
test_that("'infer_trend_cyc_seas_err_forecast' works", {
set.seed(0)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * time + age
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(time ~ RW2()) |>
set_prior(sex:time ~ Lin()) |>
fit()
mod <- set_n_draw(mod, 5)
mod <- fit(mod)
comp_est <- components(mod)
comp_forecast <- forecast(mod, labels = 2006:2007, output = "components")
dimnames_terms_forecast <- make_dimnames_terms_forecast(dimnames_terms = mod$dimnames_terms,
var_time = mod$var_time,
labels_forecast = 2006:2007,
time_only = TRUE)
ans_obtained <- infer_trend_cyc_seas_err_forecast(components = comp_forecast,
priors = mod$priors,
dimnames_terms = dimnames_terms_forecast,
var_time = mod$var_time,
var_age = mod$var_age)
expect_equal(ans_obtained[1:3], comp_forecast[1:3])
})
## 'infer_trend_cyc_seas_err_seasfix_forecast' --------------------------------
test_that("'infer_trend_cyc_seas_err_seasfix_forecast' works", {
set.seed(0)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * time + age
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ RW_Seas(n = 3)) |>
set_n_draw(n = 10) |>
fit()
components <- components(mod)
ans <- infer_trend_cyc_seas_err_seasfix_forecast(prior = mod$priors[["sex:time"]],
dimnames_term = mod$dimnames_terms[["sex:time"]],
var_time = mod$var_time,
var_age = mod$var_age,
components = components)
season <- ans$.fitted[ans$term == "sex:time" & ans$component == "season"]
trend <- ans$.fitted[ans$term == "sex:time" & ans$component == "trend"]
effect <- ans$.fitted[ans$term == "sex:time" & ans$component == "effect"]
expect_equal(effect, season + trend)
})
## 'infer_trend_cyc_seas_err_seasvary_forecast' -------------------------------
test_that("'infer_trend_cyc_seas_err_seasvary_forecast' works", {
set.seed(0)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * time + age
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ RW_Seas(n = 3, s = 1)) |>
set_n_draw(n = 10) |>
fit()
components <- components(mod)
ans <- infer_trend_cyc_seas_err_seasvary_forecast(prior = mod$priors[["sex:time"]],
dimnames_term = mod$dimnames_terms[["sex:time"]],
var_time = mod$var_time,
var_age = mod$var_age,
components = components)
season <- ans$.fitted[ans$term == "sex:time" & ans$component == "season"]
trend <- ans$.fitted[ans$term == "sex:time" & ans$component == "trend"]
effect <- ans$.fitted[ans$term == "sex:time" & ans$component == "effect"]
expect_equal(effect, season + trend)
})
## 'rvec_to_mean' -------------------------------------------------------------
test_that("'rvec_to_mean' works with valid inputs", {
data <- tibble::tibble(a = 1:3,
b = rvec::rvec(matrix(1:12, nr = 3)),
c = rvec::rvec(matrix(FALSE, nr = 3, nc = 2)),
d = c("a", "b", "c"))
ans_obtained <- rvec_to_mean(data)
ans_expected <- tibble::tibble(a = 1:3,
b = rowMeans(matrix(1:12, nr = 3)),
c = rowMeans(matrix(FALSE, nr = 3, nc = 2)),
d = c("a", "b", "c"))
expect_identical(ans_obtained, ans_expected)
})
## 'rmvnorm_chol', 'rmvnorm_eigen' --------------------------------------------
test_that("'rmvnorm_chol' and 'rmvnorm_eigen' give the same answer", {
set.seed(0)
prec <- crossprod(matrix(rnorm(25), 5))
mean <- rnorm(5)
R_prec <- chol(prec)
scaled_eigen <- make_scaled_eigen(prec)
ans_chol <- rmvnorm_chol(n = 100000, mean = mean, R_prec = R_prec)
ans_eigen <- rmvnorm_eigen(n = 100000, mean = mean, scaled_eigen = scaled_eigen)
expect_equal(rowMeans(ans_chol), rowMeans(ans_eigen), tolerance = 0.02)
expect_equal(cov(t(ans_chol)), cov(t(ans_eigen)), tolerance = 0.02)
})
## 'sort_components' ----------------------------------------------------------
test_that("'sort_components' works with valid inputs", {
components <- tibble::tribble(~term, ~component, ~level,
"sex", "effect", "m",
"time", "season", "2000",
"sex", "hyper", "sd",
"sex", "effect", "f",
"(Intercept)", "effect", "(Intercept)",
"time", "effect", "2000")
mod <- list(formula = deaths ~ time + sex)
ans_obtained <- sort_components(components = components, mod = mod)
ans_expected <- tibble::tribble(~term, ~component, ~level,
"(Intercept)", "effect", "(Intercept)",
"time", "effect", "2000",
"time", "season", "2000",
"sex", "effect", "m",
"sex", "effect", "f",
"sex", "hyper", "sd")
expect_identical(ans_obtained, ans_expected)
})
test_that("'sort_components' raises correct effor with invalid component", {
components <- tibble::tribble(~term, ~component, ~level,
"(Intercept)", "effect", "(Intercept)",
"time", "season", "2000",
"sex", "wrong", "sd")
expect_error(sort_components(components, mod = list(formula = deaths ~ age)),
"Internal error: \"wrong\" not a valid value for `component`.")
})
## 'transform_hyper_ar' -------------------------------------------------------
test_that("'transform_hyper_ar' works with 'bage_prior_ar - AR1'", {
shifted_invlogit <- function(x) {
ans <- exp(x) / (1 + exp(x))
0.18 * ans + 0.8
}
l <- transform_hyper_ar(prior = AR1())
expect_equal(l[[1]](0.35), shifted_invlogit(0.35))
expect_equal(l[[2]](0.35), exp(0.35))
})
test_that("'transform_hyper_ar' works with 'bage_prior_svd_ar - AR'", {
shifted_invlogit <- function(x) {
ans <- exp(x) / (1 + exp(x))
2 * ans - 1
}
l <- transform_hyper_ar(prior = SVD_AR(HMD))
expect_equal(l[[1]](0.35), shifted_invlogit(0.35))
expect_equal(l[[2]](0.35), shifted_invlogit(0.35))
expect_equal(l[[3]](0.35), exp(0.35))
})
Try the bage package in your browser
Any scripts or data that you put into this service are public.
bage documentation built on April 3, 2025, 8:53 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
## 'combine_stored_draws_point_inner_outer' -----------------------------------------
test_that("'combine_stored_draws_point_inner_outer' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:4,
sex = c("F", "M"),
region = c("a", "b"),
time = 2001:2005)
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + region * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, region:time ~ Lin())
mod <- set_n_draw(mod, n_draw = 10)
vars_inner <- c("age", "sex")
use_term <- make_use_term(mod = mod, vars_inner = vars_inner)
mod_inner <- reduce_model_terms(mod = mod, use_term = use_term)
mod_inner <- fit_default(mod_inner, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
mod_outer <- reduce_model_terms(mod = mod, use_term = !use_term)
mod_outer <- fit_default(mod_outer, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
mod_comb <- combine_stored_draws_point_inner_outer(mod = mod,
mod_inner = mod_inner,
mod_outer = mod_outer,
use_term = use_term)
terms <- make_terms_effectfree(mod)
is_inner <- terms %in% c("(Intercept)", "age", "sex", "age:sex")
is_outer <- terms %in% c("region", "time", "region:time")
expect_identical(mod_comb$draws_effectfree[is_inner, ],
mod_inner$draws_effectfree)
expect_identical(mod_comb$draws_effectfree[is_outer, ],
mod_outer$draws_effectfree)
expect_identical(ncol(mod_comb$draws_hyper),
ncol(mod_outer$draws_hyper))
expect_identical(ncol(mod_comb$draws_hyperrand),
ncol(mod_outer$draws_hyperrand))
expect_identical(mod_comb$point_effectfree[is_inner],
mod_inner$point_effectfree)
expect_identical(mod_comb$point_effectfree[is_outer],
mod_outer$point_effectfree)
})
## 'con_by_fitted' ------------------------------------------------------------
test_that("'con_by_fitted' works", {
set.seed(0)
prior <- RW()
fitted <- rvec::rnorm_rvec(n = 100, n_draw = 10)
along <- "time"
dimnames_term <- list(time = 2001:2010,
age = 0:4,
sex = 1:2)
var_time <- "time"
var_age <- "age"
ans <- con_by_fitted(prior = prior,
fitted = fitted,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age)
expect_equal(sum(ans[c(3, 13, 23, 33, 43)]),
sum(ans[c(4, 14, 24, 34, 44)]))
expect_equal(ans[10], -ans[60])
})
## 'draw_vals_components_fitted' ----------------------------------------------
test_that("'draw_vals_components_fitted' works", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp())
mod <- set_prior(mod, age:sex ~ SVD(HMD))
mod <- fit(mod)
set.seed(0)
ans <- draw_vals_components_fitted(mod)
expect_identical(names(ans), c("term", "component", "level", ".fitted"))
i <- ans$term == "age:sex" & ans$component == "effect"
expect_true(mean(abs(as.numeric(sum(ans$.fitted[i])))) > 0)
})
## 'generate_prior_helper' ----------------------------------------------------
test_that("'generate_prior_helper' works with valid inputs, not along-by", {
x <- N()
ans_obtained <- generate_prior_helper(x = x, n_element = 5, n_draw = 2)
draw <- factor(rep(c("Draw 1", "Draw 2"), each = 5))
ans_expected <- list(ans = tibble::tibble(draw = draw,
element = rep(1:5, times = 2)),
matrix_along_by = matrix(0:4, nc = 1, dimnames = list(seq_len(5), NULL)),
levels_effect = 1:5)
expect_identical(ans_obtained, ans_expected)
})
test_that("'generate_prior_helper' works with valid inputs, along-by, n_by = 1", {
x <- AR()
ans_obtained <- generate_prior_helper(x = x, n_along = 5, n_by = 1, n_draw = 2)
draw <- factor(rep(c("Draw 1", "Draw 2"), each = 5))
by <- factor(rep("By 1", 10))
levels_effect <- paste(1, 1:5, sep = ".")
ans_expected <- list(ans = tibble::tibble(draw = draw,
by = by,
along = rep(1:5, times = 2)),
matrix_along_by = matrix(0:4, nc = 1, dimnames = list(seq_len(5), 1L)),
levels_effect = levels_effect)
expect_identical(ans_obtained, ans_expected)
})
test_that("'generate_prior_helper' works with valid inputs, along-by, n_by = 2", {
x <- AR()
ans_obtained <- generate_prior_helper(x = x, n_along = 5, n_by = 2, n_draw = 2)
draw <- factor(rep(c("Draw 1", "Draw 2"), each = 10))
by <- factor(rep(rep(paste("By", 1:2), each = 5), times = 2))
ans_expected <- list(ans = tibble::tibble(draw = draw,
by = by,
along = rep(1:5, times = 4)),
matrix_along_by = matrix(0:9, nc = 2, dimnames = list(1:5, 1:2)),
levels_effect = paste(rep(1:2, each = 5), 1:5, sep = "."))
expect_identical(ans_obtained, ans_expected)
})
## 'generate_prior_svd_helper' ------------------------------------------------
test_that("'generate_prior_svd_helper' works with valid inputs - n_element", {
x <- SVD(LFP)
set.seed(0)
n_element <- 20
n_draw <- 25
ans_obtained <- generate_prior_svd_helper(x, n_element = n_element, n_draw = n_draw)
expect_true(is.list(ans_obtained))
})
test_that("'generate_prior_svd_helper' works with valid inputs - n_by = 1", {
x <- SVD_RW(LFP)
set.seed(0)
n_along <- 20
n_by <- 2
n_draw <- 25
ans_obtained <- generate_prior_svd_helper(x, n_along = n_along, n_by = n_by, n_draw = n_draw)
expect_true(is.list(ans_obtained))
})
## 'generate_ssvd_helper' -----------------------------------------------------------------
test_that("'generate_ssvd_helper' works with valid inputs - indep = NULL, n_element = 1", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_element = 1,
n_draw = 3,
n_comp = 2,
indep = NULL,
age_labels = NULL)
expect_identical(nrow(ans$matrix), length(unique(ans$ans$age)))
expect_identical(ncol(ans$matrix), 2L)
})
test_that("'generate_ssvd_helper' works with valid inputs - indep = TRUE, n_by = 1, n_along = 3", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_along = 3,
n_by = 1,
n_draw = 3,
n_comp = 2,
indep = TRUE,
age_labels = NULL)
expect_identical(nrow(ans$matrix), 3L * nrow(unique(ans$ans[c("age", "sexgender")])))
expect_identical(ncol(ans$matrix), 12L)
})
test_that("'generate_ssvd_helper' works with valid inputs - indep = FALSE, n_by = 1", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_along = 2,
n_by = 1,
n_draw = 3,
n_comp = 2,
indep = FALSE,
age_labels = NULL)
expect_identical(nrow(ans$matrix), 2L * nrow(unique(ans$ans[c("age", "sexgender")])))
expect_identical(ncol(ans$matrix), 4L)
})
test_that("'generate_ssvd_helper' works with valid inputs - indep = TRUE, n_by = 2", {
set.seed(0)
ans <- generate_ssvd_helper(ssvd = LFP,
n_along = 3,
n_by = 2,
n_draw = 3,
n_comp = 2,
indep = TRUE,
age_labels = NULL)
expect_identical(nrow(ans$matrix), 3L * nrow(unique(ans$ans[c("by", "age", "sexgender")])))
expect_identical(ncol(ans$matrix), 24L)
expect_identical(nrow(ans$matrix_along_by), 3L)
})
## 'get_term_from_est' --------------------------------------------------------
test_that("'get_term_from_est' works", {
est <- list(effectfree = c(a = 1, a = 2, b = 3, c = 4, c = 5),
hyper = c(a = 1),
hyperrandfree = double(),
log_disp = c(disp = 3))
ans_obtained <- get_term_from_est(est = est, index_term = c(1L, 4L, 7L))
ans_expected <- c("a", "c", "disp")
expect_identical(ans_obtained, ans_expected)
})
## 'get_disp' -----------------------------------------------------------------
test_that("'get_disp' works - unfitted", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 5)
set.seed(1)
ans_obtained <- get_disp(mod)
set.seed(1)
ans_expected <- draw_vals_disp(mod, n_sim = 5)
expect_identical(ans_obtained, ans_expected)
})
test_that("'get_disp' works - fitted", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 5)
mod <- fit(mod)
set.seed(1)
ans_obtained <- get_disp(mod)
set.seed(1)
ans_expected <- rvec::rvec(matrix(mod$draws_disp, nr = 1))
expect_identical(ans_obtained, ans_expected)
})
## 'insert_after' -------------------------------------------------------
test_that("'insert_after' works with data frames", {
df <- data.frame(x = 1:3, y = 3:1)
x <- 11:13
nm_x = "new"
ans_obtained <- insert_after(df = df,
nm_after = "x",
x = x,
nm_x = nm_x)
ans_expected <- data.frame(x = 1:3, new = 11:13, y = 3:1)
expect_identical(ans_obtained, ans_expected)
ans_obtained <- insert_after(df = df,
nm_after = "y",
x = x,
nm_x = nm_x)
ans_expected <- data.frame(x = 1:3, y = 3:1, new = 11:13)
expect_identical(ans_obtained, ans_expected)
})
test_that("'insert_after' works with tibbles", {
df <- tibble::tibble(x = 1:3, y = 3:1)
x <- 11:13
nm_x = "new"
ans_obtained <- insert_after(df = df,
nm_after = "x",
x = x,
nm_x = nm_x)
ans_expected <- tibble(x = 1:3, new = 11:13, y = 3:1)
expect_identical(ans_obtained, ans_expected)
ans_obtained <- insert_after(df = df,
nm_after = "y",
x = x,
nm_x = nm_x)
ans_expected <- tibble(x = 1:3, y = 3:1, new = 11:13)
expect_identical(ans_obtained, ans_expected)
})
## 'is_same_class' ------------------------------------------------------------
test_that("'is_same_class' returns TRUE when classes same", {
expect_true(is_same_class(AR1(), AR1()))
expect_true(is_same_class(1L, 2L))
})
test_that("'is_same_class' returns FALSE when classes different", {
expect_false(is_same_class(AR1(), N()))
expect_false(is_same_class(1L, FALSE))
})
## 'make_combined_matrix_effect_outcome' -----------------------------------------
test_that("'make_combined_matrix_effect_outcome' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_combined_matrix_effect_outcome(mod)
expect_identical(nrow(ans_obtained), nrow(data))
expect_identical(ncol(ans_obtained), length(make_terms_effects(mod$dimnames_terms)))
expect_false(all(ans_obtained[,1] == 0))
})
## 'make_combined_matrix_effectfree_effect' -----------------------------------------
test_that("'make_combined_matrix_effectfree_effect' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_combined_matrix_effectfree_effect(mod)
expect_identical(nrow(ans_obtained), length(make_terms_effects(mod$dimnames_terms)))
expect_identical(ncol(ans_obtained), length(make_terms_effectfree(mod)))
})
## 'make_combined_offset_effectfree_effect' -----------------------------------------
test_that("'make_combined_offset_effectfree_effect' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_combined_offset_effectfree_effect(mod)
expect_identical(length(ans_obtained), length(make_terms_effects(mod$dimnames_terms)))
expect_true(all(ans_obtained == 0))
})
## 'make_comp_components' -----------------------------------------------------
test_that("'make_comp_components' works - no hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
draws <- make_draws_components(mod)
ans <- make_comp_components(mod)
expect_identical(length(draws), length(ans))
expect_setequal(ans, c("effect", "hyper", "disp"))
})
test_that("'make_comp_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ Lin())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
draws <- make_draws_components(mod)
ans <- make_comp_components(mod)
expect_identical(length(draws), length(ans))
expect_setequal(ans, c("effect", "hyper", "trend", "error", "disp"))
})
test_that("'make_comp_components' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD))
mod <- set_prior(mod, age:sex ~ Sp())
set.seed(0)
mod <- fit(mod)
term <- make_term_components(mod)
ans <- make_comp_components(mod)
expect_identical(length(term), length(ans))
expect_setequal(ans, c("effect", "hyper", "svd", "spline", "disp"))
})
## 'make_comp_hyperrand' ------------------------------------------------------
test_that("'make_comp_hyperrand' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex:time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
ans_obtained <- make_comp_hyperrand(mod)
ans_expected <- c("hyper", "hyper", rep(c("trend", "error"), each = 12))
expect_identical(ans_obtained, ans_expected)
})
## 'make_copies_repdata' ------------------------------------------------------
test_that("'make_copies_repdata' works with valid inputs", {
data <- tibble::tibble(age = 0:4, deaths = 1:5, population = 11:15)
n <- 2
ans_obtained <- make_copies_repdata(data = data, n = n)
ans_expected <- tibble::tibble(.replicate = factor(rep(c("Original", "Replicate 1", "Replicate 2"),
each = 5),
levels = c("Original",
"Replicate 1", "Replicate 2")),
age = rep(0:4, 3),
deaths = rep(1:5, 3),
population = rep(11:15, 3))
expect_identical(ans_obtained, ans_expected)
})
## 'make_draws_components' -----------------------------------------------
test_that("'make_draws_components' works - no svd, spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
expect_identical(length(ans_obtained),
length(make_effectfree(mod)) +
length(make_hyper(mod)) +
length(make_hyperrand(mod)) + 1L)
})
test_that("'make_draws_components' works - has spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 6))
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
expect_identical(length(ans_obtained),
length(make_effectfree(mod)) + length(unique(data$age)) +
length(make_hyper(mod)) +
length(make_hyperrand(mod)) + 1L)
})
test_that("'make_draws_components' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD, n_comp = 3))
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
expect_identical(length(ans_obtained),
length(make_effectfree(mod)) + length(unique(data$age)) +
length(make_hyper(mod)) +
length(make_hyperrand(mod)) + 1L)
})
test_that("'make_draws_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age:sex ~ Lin())
mod <- fit(mod)
set.seed(0)
ans_obtained <- make_draws_components(mod)
expect_true(rvec::is_rvec(ans_obtained))
})
## 'make_draws_disp' ----------------------------------------------------
test_that("'make_draws_disp' works", {
set.seed(0)
draws_post <- matrix(rnorm(13 * 5), nrow = 13)
ans_obtained <- make_draws_disp(draws_post)
ans_expected <- exp(draws_post[13,])
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_draws_effectfree' ----------------------------------------------------
test_that("'make_draws_effectfree' works", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = rnorm(2),
hyperrand = numeric(),
log_disp = 0.5)
draws_post <- matrix(rnorm(13 * 5), nrow = 13)
ans_obtained <- make_draws_effectfree(est = est, draws_post = draws_post)
ans_expected <- draws_post[1:10,]
expect_identical(ans_obtained, ans_expected)
})
## 'make_draws_hyper' ----------------------------------------------------
test_that("'make_draws_hyper' works", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = runif(4),
hyperrand = numeric(),
log_disp = 0.5)
transforms_hyper <- list(identity, identity, exp, exp)
draws_post <- matrix(rnorm(15 * 5), nrow = 15)
ans_obtained <- make_draws_hyper(est = est,
transforms_hyper = transforms_hyper,
draws_post = draws_post)
ans_expected <- draws_post[11:14, ]
ans_expected[3:4,] <- exp(ans_expected[3:4,])
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_draws_hyperrandfree' -------------------------------------------------
test_that("'make_draws_hyperrandfree' works - has hyperrandfree", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = runif(4),
hyperrandfree = rnorm(5),
log_disp = 0.5)
draws_post <- matrix(rnorm(20 * 5), nrow = 20)
ans_obtained <- make_draws_hyperrandfree(est = est,
draws_post = draws_post)
ans_expected <- draws_post[15:19, ]
expect_identical(unname(ans_obtained), ans_expected)
})
test_that("'make_draws_hyperrandfree' works - no hyperrandfree", {
set.seed(0)
est <- list(effectfree = rnorm(10),
hyper = runif(4),
hyperrandfree = numeric(),
log_disp = 0.5)
draws_post <- matrix(rnorm(15 * 5), nrow = 15)
ans_obtained <- make_draws_hyperrandfree(est = est,
draws_post = draws_post)
ans_expected <- matrix(0, nrow = 0, ncol = 5)
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_draws_post' ------------------------------------------------------
test_that("'make_draws_post' works with valid inputs - has R_prec", {
set.seed(0)
est <- list(effectfree = c("(Intercept)" = -3,
sex = c(-1, 1),
age = rnorm(10)),
log_disp = 0.2)
prec <- diag(12)
prec <- Matrix::Matrix(prec)
map <- list(effectfree = factor(c(1, NA, NA, 2:11)),
disp = factor(1))
ans <- make_draws_post(est = est,
prec = prec,
map = map,
n_draw = 500000L)
expect_equal(rowMeans(ans), unlist(est, use.names = FALSE), tolerance = 0.02)
expect_equal(solve(cov(t(ans[c(1, 4:14),]))), as.matrix(prec), tolerance = 0.02)
})
## 'fit_default' --------------------------------------------------------------
## lots more tests for 'fit' method
test_that("'fit_default' works with pois, optimzier is 'multi'", {
set.seed(10)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(age ~ AR()) |>
set_prior(age:sex ~ RW2(sd = 0)) |>
set_prior(age:sex:time ~ AR())
ans_obtained <- fit_default(mod, optimizer = "multi", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
expect_s3_class(ans_obtained, "bage_mod")
})
test_that("'fit_default' works with pois, optimzier is 'nlminb'", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
expect_s3_class(ans_obtained, "bage_mod")
})
test_that("'fit_default' works with pois - start_oldpar", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = FALSE)
ans_obtained <- fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = TRUE)
expect_s3_class(ans_obtained, "bage_mod")
})
test_that("'fit_default' gives error with 'start_oldpar' if model fitted", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
expect_error(fit_default(mod, optimizer = "nlminb", quiet = TRUE, aggregate = TRUE,
start_oldpar = TRUE),
"`start_oldpar` is TRUE but model has not been fitted.")
})
## 'fit_inner_outer' ----------------------------------------------------------
test_that("'fit_inner_outer' works with with pois", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$popn <- rpois(n = nrow(data), lambda = 1000)
data$deaths <- rnbinom(n = nrow(data), mu = 0.1 * data$popn, size = 100)
formula <- deaths ~ age * sex + region * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
set.seed(0)
ans_inner_outer <- fit_inner_outer(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
vars_inner = c("age", "sex"))
set.seed(0)
ans_default <- fit_default(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
aggregate = TRUE)
aug_inner_outer <- ans_inner_outer |>
augment()
fit_inner_outer <- rvec::draws_median(aug_inner_outer$.fitted)
aug_default <- ans_default |>
augment()
fit_default <- rvec::draws_median(aug_default$.fitted)
expect_true(cor(fit_inner_outer, fit_default) > 0.98)
})
test_that("'fit_inner_outer' works with with binom", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$popn <- rpois(n = nrow(data), lambda = 1000)
data$deaths <- rbinom(n = nrow(data), prob = 0.2, size = data$popn)
formula <- deaths ~ age * sex + region * time
mod <- mod_binom(formula = formula,
data = data,
size = popn)
set.seed(0)
ans_inner_outer <- fit_inner_outer(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
vars_inner = NULL)
set.seed(0)
ans_default <- fit_default(mod,
optimizer = "nlminb",
quiet = TRUE,
start_oldpar = FALSE,
aggregate = TRUE)
aug_inner_outer <- ans_inner_outer |>
augment()
fit_inner_outer <- rvec::draws_median(aug_inner_outer$.fitted)
aug_default <- ans_default |>
augment()
fit_default <- rvec::draws_median(aug_default$.fitted)
expect_true(cor(fit_inner_outer, fit_default) > 0.98)
})
test_that("'fit_inner_outer' works with with norm", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$wt <- rpois(n = nrow(data), lambda = 10)
data$income <- rnorm(n = nrow(data), mean = data$age + data$time/100, sd = 5 / sqrt(data$wt))
formula <- income ~ age * sex + region * time
mod <- mod_norm(formula = formula,
data = data,
weights = wt)
set.seed(0)
ans_inner_outer <- fit_inner_outer(mod,
optimizer = "BFGS",
quiet = TRUE,
start_oldpar = FALSE,
vars_inner = c("age", "sex"))
set.seed(0)
ans_default <- fit_default(mod,
optimizer = "BFGS",
quiet = TRUE,
start_oldpar = FALSE,
aggregate = TRUE)
aug_inner_outer <- ans_inner_outer |>
augment()
fit_inner_outer <- rvec::draws_median(aug_inner_outer$.fitted)
aug_default <- ans_default |>
augment()
fit_default <- rvec::draws_median(aug_default$.fitted)
expect_true(cor(fit_inner_outer, fit_default) > 0.98)
})
test_that("'fit_inner_outer' throws error when 'start_oldpar' is TRUE", {
set.seed(0)
data <- expand.grid(age = 0:5,
time = 2000:2003,
sex = c("F", "M"),
region = c("a", "b"))
data$wt <- rpois(n = nrow(data), lambda = 10)
data$income <- rnorm(n = nrow(data), mean = data$age + data$time/100, sd = 5 / sqrt(data$wt))
formula <- income ~ age * sex + region * time
mod <- mod_norm(formula = formula,
data = data,
weights = wt)
set.seed(0)
expect_error(fit_inner_outer(mod,
optimizer = "BFGS",
quiet = TRUE,
start_oldpar = TRUE,
vars_inner = c("age", "sex")),
"`start_oldpar` must be FALSE when using \"inner-outer\" method.")
})
## 'make_along_mod' -----------------------------------------------------------
test_that("'make_along_mod' works", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:time ~ Sp(n_comp = 5, along = "age"))
ans_obtained <- make_along_mod(mod)
ans_expected <- c("(Intercept)" = NA,
sex = NA,
age = "age",
time = "time",
"sex:age" = "age",
"age:time" = "age")
expect_identical(ans_obtained, ans_expected)
})
## 'make_hyperrand' -----------------------------------------------------------
test_that("'make_hyperrand' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ Lin())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
ans <- make_hyperrand(mod)
expect_identical(length(ans), 26L)
})
test_that("'make_hyperrand' works - no hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ AR())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
ans <- make_hyperrand(mod)
expect_identical(length(ans), 0L)
})
## 'make_hyperrand_lin' -------------------------------------------------------
test_that("'make_hyperrand_lin' works with main effect", {
prior <- Lin()
hyperrandfree <- rvec::rnorm_rvec(n = 1, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_lin(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
trend <- hyperrandfree * (1:10 - mean(1:10))
error <- effectfree - trend
slope <- trend[2] - trend[1]
ans_expected <- c(slope, trend, error)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_lin' works with interaction, con is 'none'", {
prior <- Lin_AR()
hyperrandfree <- rvec::rnorm_rvec(n = 2, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_lin(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
trend <- c(hyperrandfree[1] * (1:10 - mean(1:10)),
hyperrandfree[2] * (1:10 - mean(1:10)))
error <- effectfree - trend
slope <- c(trend[2] - trend[1], trend[12] - trend[11])
ans_expected <- c(slope, trend, error)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_lin' works with interaction, con is 'by'", {
prior <- Lin_AR(con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 1, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_lin(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
trend <- c(-sqrt(0.5) * hyperrandfree * (1:10 - mean(1:10)),
sqrt(0.5) * hyperrandfree * (1:10 - mean(1:10)))
error <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - trend
slope <- c(trend[2] - trend[1], trend[12] - trend[11])
ans_expected <- c(slope, trend, error)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_randomseasfix' -----------------------------------------------
test_that("'make_hyperrand_randomseasfix' works with main effect", {
prior <- RW_Seas(n = 4)
hyperrandfree <- rvec::rnorm_rvec(n = 3, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(hyperrandfree[1],
hyperrandfree[2],
hyperrandfree[3],
-sum(hyperrandfree))[c(1:4, 1:4, 1:2)]
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasfix' works with interaction, random_sum is FALSE", {
set.seed(0)
prior <- RW2_Seas(n = 4)
hyperrandfree <- rvec::rnorm_rvec(n = 6, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(hyperrandfree[1],
hyperrandfree[2],
hyperrandfree[3],
-sum(hyperrandfree[1:3]))[c(1:4, 1:4, 1:2)]
season2 <- c(hyperrandfree[4],
hyperrandfree[5],
hyperrandfree[6],
-sum(hyperrandfree[4:6]))[c(1:4, 1:4, 1:2)]
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasfix' works with interaction, con is 'by'", {
prior <- RW_Seas(n = 4, con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 3, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(hyperrandfree[1],
hyperrandfree[2],
hyperrandfree[3],
-sum(hyperrandfree[1:3]))[c(1:4, 1:4, 1:2)]
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_randomseasvary' ----------------------------------------------
test_that("'make_hyperrand_randomseasvary' works with main effect", {
prior <- RW_Seas(n = 4, s_seas = 1)
hyperrandfree <- rvec::rnorm_rvec(n = 8, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(hyperrandfree[1:3],
-sum(hyperrandfree[1:3]),
hyperrandfree[4:6],
-sum(hyperrandfree[4:6]),
hyperrandfree[7:8])
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasvary' works with interaction, con is 'none'", {
set.seed(0)
prior <- RW2_Seas(n = 4, s_seas = 1)
hyperrandfree <- rvec::rnorm_rvec(n = 16, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(hyperrandfree[1:3],
-sum(hyperrandfree[1:3]),
hyperrandfree[4:6],
-sum(hyperrandfree[4:6]),
hyperrandfree[7:8])
season2 <- c(hyperrandfree[9:11],
-sum(hyperrandfree[9:11]),
hyperrandfree[12:14],
-sum(hyperrandfree[12:14]),
hyperrandfree[15:16])
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_randomseasvary' works with interaction, con is 'by'", {
prior <- RW_Seas(n_seas = 4, s_seas = 1, con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 8, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_randomseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(hyperrandfree[1:3],
-sum(hyperrandfree[1:3]),
hyperrandfree[4:6],
-sum(hyperrandfree[4:6]),
hyperrandfree[7:8])
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_zeroseasfix' -----------------------------------------------
test_that("'make_hyperrand_zeroseasfix' works with main effect", {
prior <- RW_Seas(n = 4, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 2, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree))[c(1:4, 1:4, 1:2)]
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasfix' works with interaction, con is 'none'", {
set.seed(0)
prior <- RW2_Seas(n = 4, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 4, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]))[c(1:4, 1:4, 1:2)]
season2 <- c(effectfree[11],
effectfree[11] + hyperrandfree[3],
effectfree[11] + hyperrandfree[4],
-3 * effectfree[11] - sum(hyperrandfree[3:4]))[c(1:4, 1:4, 1:2)]
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasfix' works with interaction, con is 'by'", {
prior <- RW_Seas(n = 4, sd = 0, con = "by")
hyperrandfree <- rvec::rnorm_rvec(n = 2, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasfix(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]))[c(1:4, 1:4, 1:2)]
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_hyperrand_zeroseasvary' ----------------------------------------------
test_that("'make_hyperrand_zeroseasvary' works with main effect", {
prior <- RW_Seas(n = 4, s_seas = 1, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 7, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10, n_draw = 10)
dimnames_term <- list(time = 2001:2010)
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]),
effectfree[1] + hyperrandfree[3:5],
-3 * effectfree[1] - sum(hyperrandfree[3:5]),
effectfree[1] + hyperrandfree[6:7])
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasvary' works with interaction, con is 'none'", {
set.seed(0)
prior <- RW2_Seas(n = 4, s_seas = 1, sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 14, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 20, n_draw = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season1 <- c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]),
effectfree[1] + hyperrandfree[3:5],
-3 * effectfree[1] - sum(hyperrandfree[3:5]),
effectfree[1] + hyperrandfree[6:7])
season2 <- c(effectfree[11],
effectfree[11] + hyperrandfree[8],
effectfree[11] + hyperrandfree[9],
-3 * effectfree[11] - sum(hyperrandfree[8:9]),
effectfree[11] + hyperrandfree[10:12],
-3 * effectfree[11] - sum(hyperrandfree[10:12]),
effectfree[11] + hyperrandfree[13:14])
season <- c(season1, season2)
trend <- effectfree - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_hyperrand_zeroseasvary' works with interaction, con is 'by'", {
prior <- RW_Seas(n_seas = 4, s_seas = 1, con = "by", sd = 0)
hyperrandfree <- rvec::rnorm_rvec(n = 7, n_draw = 10)
effectfree <- rvec::rnorm_rvec(n = 10)
dimnames_term <- list(time = 2001:2010, sex = c("f", "m"))
var_time <- "time"
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_hyperrand_zeroseasvary(prior = prior,
hyperrandfree = hyperrandfree,
effectfree = effectfree,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
season <- vctrs::vec_c(effectfree[1],
effectfree[1] + hyperrandfree[1],
effectfree[1] + hyperrandfree[2],
-3 * effectfree[1] - sum(hyperrandfree[1:2]),
effectfree[1] + hyperrandfree[3:5],
-3 * effectfree[1] - sum(hyperrandfree[3:5]),
effectfree[1] + hyperrandfree[6:7])
season <- c(-sqrt(0.5) * season,
sqrt(0.5) * season)
trend <- c(-sqrt(0.5) * effectfree, sqrt(0.5) * effectfree) - season
ans_expected <- c(trend, season)
expect_equal(ans_obtained, ans_expected)
})
## 'make_levels_spline' -------------------------------------------------------
test_that("'make_levels_spline' works - unlist is FALSE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:sex ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:time ~ Sp(n_comp = 5, along = "age"))
set.seed(0)
ans_obtained <- make_levels_spline(mod, unlist = FALSE)
ans_expected <- list("(Intercept)" = NULL,
sex = NULL,
age = paste0("comp", 1:5),
time = NULL,
"sex:age" = paste(c("F", "M"), rep(paste0("comp", 1:5), each = 2), sep = "."),
"age:time" = paste0("comp", 1:5, ".", rep(2000:2005, each = 5)))
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_spline' works - unlist is TRUE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:sex ~ Sp(n_comp = 5))
mod <- set_prior(mod, age:time ~ Sp(n_comp = 5, along = "age"))
set.seed(0)
ans_obtained <- make_levels_spline(mod, unlist = TRUE)
ans_expected <- c(paste0("comp", 1:5),
paste(c("F", "M"), rep(paste0("comp", 1:5), each = 2), sep = "."),
paste(paste0("comp", 1:5), rep(2000:2005, each = 5), sep = "."))
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_spline_term' --------------------------------------------------
test_that("'make_levels_spline_term' works - con is 'none'", {
prior <- Sp(n = 5)
dimnames_term <- list(reg = 1:2,
age = 1:20)
ans_obtained <- make_levels_spline_term(prior = prior,
dimnames_term = dimnames_term,
var_time = "time",
var_age = "age")
ans_expected <- paste(1:2, rep(paste0("comp", 1:5), each = 2), sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_spline_term' works - con is 'by'", {
prior <- Sp(n = 5, con = "by")
dimnames_term <- list(reg = 1:2,
age = 1:20)
ans_obtained <- make_levels_spline_term(prior = prior,
dimnames_term = dimnames_term,
var_time = "time",
var_age = "age")
ans_expected <- paste("reg1", paste0("comp", 1:5), sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_spline_term' works - con is 'by', more dimensions", {
prior <- Sp(n = 5, con = "by")
dimnames_term <- list(reg = 1:4,
age = 1:20,
sex = c("f", "m"))
ans_obtained <- make_levels_spline_term(prior = prior,
dimnames_term = dimnames_term,
var_time = "time",
var_age = "age")
ans_expected <- paste(paste0("reg", 1:3),
rep(paste0("comp", 1:5), each = 3),
rep("sex1", times = 15),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_svd' ----------------------------------------------------------
test_that("'make_levels_svd' works - unlist is FALSE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD))
mod <- set_prior(mod, age:sex ~ SVD(HMD))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD))
set.seed(0)
ans_obtained <- make_levels_svd(mod, unlist = FALSE)
ans_expected <- list("(Intercept)" = NULL,
sex = NULL,
age = paste0("comp", 1:3),
time = NULL,
"sex:age" = paste0(rep(c("F", "M"), each = 3), ".comp", 1:3),
"age:time" = paste0("comp", 1:3, ".", rep(2000:2005, each = 3)))
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd' works - unlist is TRUE", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * age + age*time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:sex ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD, n_comp = 5))
set.seed(0)
ans_obtained <- make_levels_svd(mod, unlist = TRUE)
ans_expected <- c(paste0("comp", 1:5),
paste0(rep(c("F", "M"), each = 5), ".comp", 1:5),
paste0("comp", 1:5, ".", rep(2000:2005, each = 5)))
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_svd_term' ----------------------------------------------------------
test_that("'make_levels_svd_term' works - total, main effect", {
prior <- SVD(HMD)
dimnames_term <- list(age = c(0:59, "60+"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste0("comp", 1:3)
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - joint, age:sex", {
prior <- SVD(HMD, n_comp = 5)
dimnames_term <- list(age = c(0:59, "60+"), sex = c("M", "F"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(rep(c("M", "F"), each = 5), paste0("comp", 1:5), sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - indep, age:sex:reg", {
prior <- SVD(HMD, indep = FALSE, n_comp = 5)
dimnames_term <- list(reg = 1:2, age = c(0:59, "60+"), sex = c("M", "F"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(paste0("comp", 1:5),
rep(1:2, each = 5),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - indep, age:sex:time, con is 'none'", {
prior <- SVD_RW(HMD, indep = FALSE, n_comp = 5)
dimnames_term <- list(time = 2001:2010, age = c(0:59, "60+"), sex = c("M", "F"))
var_age <- "age"
var_sexgender <- "sex"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(paste0("comp", 1:5),
rep(2001:2010, each = 5),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_levels_svd_term' works - indep, age:time:reg, con is 'by'", {
prior <- SVD_RW(HMD, indep = FALSE, n_comp = 5, con = "by")
dimnames_term <- list(time = 2001:2010, age = c(0:59, "60+"), reg = c("a", "b", "c"))
var_age <- "age"
var_sexgender <- "sex"
var_time <- "time"
ans_obtained <- make_levels_svd_term(prior = prior,
dimnames_term = dimnames_term,
var_time = var_time,
var_age = var_age,
var_sexgender = var_sexgender)
ans_expected <- paste(paste0("comp", 1:5),
rep(2001:2010, each = 5),
rep(c("reg1", "reg2"), each = 50),
sep = ".")
expect_identical(ans_obtained, ans_expected)
})
## 'make_lin_trend' -----------------------------------------------------------
test_that("'make_lin_trend' works with valid inputs - n_by = 1", {
slope <- rvec::rvec(matrix(as.numeric(1:5), nr = 1))
matrix_along_by <- matrix(0:9, nr = 10)
ans_obtained <- make_lin_trend(slope = slope,
matrix_along_by = matrix_along_by)
ones <- rep(1, times = 10)
s <- 1:10
ans_expected <- rvec::rvec(outer(ones, -0.5 * 11 * (1:5)) + outer(s, 1:5))
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_lin_trend' works with valid inputs - n_by = 2, transposed", {
slope <- rvec::rvec(matrix(rnorm(10), nr = 2))
intercept <- -0.5 * 6 * slope
matrix_along_by <- t(matrix(0:9, nr = 2))
ans_obtained <- make_lin_trend(slope = slope,
matrix_along_by = matrix_along_by)
s <- 1:5
ans_expected <- c(intercept[1] + slope[1] * s, intercept[2] + slope[2] * s)
ans_expected <- ans_expected[c(1, 6, 2, 7, 3, 8, 4, 9, 5, 10)]
expect_identical(ans_obtained, ans_expected)
})
## 'make_stored_draws' --------------------------------------------------------
test_that("'make_stored_draws' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex ~ Known(c(0.1, -0.1)))
mod <- set_n_draw(mod, n = 10)
est <- list(effectfree = c(rnorm(11), 0.1, -0.1),
hyper = rnorm(1),
hyperrandfree = numeric(),
disp = runif(1))
prec <- crossprod(matrix(rnorm(169), nr = 13))
map <- make_map(mod)
ans <- make_stored_draws(mod = mod,
est = est,
prec = prec,
map = map)
expect_identical(ncol(ans$draws_effectfree), 10L)
expect_identical(ncol(ans$draws_hyper), 10L)
expect_identical(length(ans$draws_disp), 10L)
ans2 <- make_stored_draws(mod = mod,
est = est,
prec = prec,
map = map)
expect_identical(ans, ans2)
})
## 'make_stored_point' --------------------------------------------------------
test_that("'make_stored_draws' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex ~ Known(c(0.1, -0.1)))
est <- list(effectfree = c(rnorm(11), 0.1, -0.1),
hyper = rnorm(1),
hyperrandfree = numeric(),
log_disp = runif(1))
ans <- make_stored_point(mod = mod,
est = est)
expect_identical(ans$point_effectfree, est$effectfree)
expect_identical(ans$point_hyper, exp(est$hyper))
expect_identical(ans$point_hyperrandfree, double())
expect_identical(ans$point_disp, exp(est$log_disp))
})
## 'make_par_disp' ------------------------------------------------------------
test_that("'make_par_disp' works with bage_mod_pois", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + time + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 10)
mod <- fit(mod)
components <- components(mod)
meanpar <- exp(make_linpred_comp(components = components,
data = mod$data,
dimnames_terms = mod$dimnames_terms))
disp <- components$.fitted[components$component == "disp"]
set.seed(1)
ans_obtained <- make_par_disp(mod,
meanpar = meanpar,
disp = disp)
set.seed(1)
ans_expected <- rvec::rgamma_rvec(n = length(meanpar),
data$deaths + 1/disp,
data$popn + 1/(disp*meanpar))
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_par_disp' works with bage_mod_binom", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rbinom(n = nrow(data), size = data$popn, prob = 0.3)
formula <- deaths ~ age + time + sex
mod <- mod_binom(formula = formula,
data = data,
size = popn)
mod <- fit(mod)
mod <- set_n_draw(mod, n = 10)
components <- components(mod)
invlogit <- function(x) 1 / (1 + exp(-x))
meanpar <- invlogit(make_linpred_comp(components = components,
data = mod$data,
dimnames_terms = mod$dimnames_terms))
disp <- components$.fitted[components$component == "disp"]
set.seed(1)
ans_obtained <- make_par_disp(mod,
meanpar = meanpar,
disp = disp)
set.seed(1)
ans_expected <- rvec::rbeta_rvec(n = length(meanpar),
data$deaths + meanpar/disp,
data$popn - data$deaths + (1 - meanpar)/disp)
expect_equal(ans_obtained, ans_expected)
})
## 'make_is_fixed' ------------------------------------------------------------
test_that("'make_is_fixed' works when nothing fixed", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- fit(mod)
est <- list(effectfree = make_effectfree(mod),
hyper = make_hyper(mod),
hyperrand = make_hyperrand(mod),
log_disp = 0)
map <- make_map(mod)
expect_true(is.null(map))
ans_obtained <- make_is_fixed(est = est, map = map)
ans_expected <- rep(FALSE, times = length(unlist(est)))
})
test_that("'make_is_fixed' works when Known prior", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex ~ Known(c(0.1, -0.1)))
mod <- fit(mod)
est <- list(effectfree = make_effectfree(mod),
hyper = make_hyper(mod),
hyperrand = make_hyperrand(mod),
log_disp = 0)
map <- make_map(mod)
ans_obtained <- make_is_fixed(est = est, map = map)
ans_expected <- rep(c(FALSE, TRUE, FALSE),
times = c(11,
2,
20 + 6 + length(est$hyper) +
+ length(est$log_disp)))
expect_identical(unname(ans_obtained), ans_expected)
})
## 'make_level_components' ----------------------------------------------------
test_that("'make_level_components' works - no hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_disp(mod, mean = 0)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_level_components(mod)
expect_identical(length(ans), length(comp))
})
test_that("'make_level_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_disp(mod, mean = 0)
mod <- set_prior(mod, sex:time ~ Lin())
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_level_components(mod)
expect_identical(length(ans), length(comp))
})
## 'make_levels_hyper' --------------------------------------------------------
test_that("'make_levels_hyper' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
ans_obtained <- make_levels_hyper(mod)
ans_expected <- c(age = "sd",
time = "sd",
"age:sex" = "sd")
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_hyperrand' ----------------------------------------------------
test_that("'make_levels_hyperrand' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex:time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
ans_obtained <- make_levels_hyperrand(mod, unlist = TRUE)
ans_expected <- c("slope.F", "slope.M",
rep(paste(c("F", "M"), rep(2000:2005, each = 2), sep = "."), 2))
expect_identical(ans_obtained, ans_expected)
ans_obtained <- make_levels_hyperrand(mod, unlist = FALSE)
ans_expected <- list(character(),
character(),
c("slope.F", "slope.M",
rep(paste(c("F", "M"), rep(2000:2005, each = 2), sep = "."), 2)))
expect_identical(ans_obtained, ans_expected)
})
## 'make_levels_replicate' ----------------------------------------------------
test_that("'make_levels_replicate' works", {
ans_obtained <- make_levels_replicate(n = 2, n_row_data = 3)
ans_expected <- factor(rep(c("Original", "Replicate 1", "Replicate 2"),
each = 3),
levels = c("Original", "Replicate 1", "Replicate 2"))
expect_identical(ans_obtained, ans_expected)
})
## 'make_linpred_comp' --------------------------------------------------------
test_that("'make_linpred_comp' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 10L)
comp <- components(mod, quiet = TRUE)
ans <- make_linpred_comp(components = comp,
data = mod$data,
dimnames_terms = mod$dimnames_terms)
expect_identical(length(ans), length(mod$outcome))
})
## 'make_linpred_raw' ---------------------------------------------------------
test_that("'make_linpred_raw' works with valid inputs - point is FALSE", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 10L)
mod <- fit(mod)
ans_obtained <- make_linpred_raw(mod, point = FALSE)
comp <- components(mod, quiet = TRUE)
ans_expected <- make_linpred_comp(components = comp,
data = mod$data,
dimnames_terms = mod$dimnames_terms)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_linpred_raw' works with valid inputs - point is TRUE", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n_draw = 10L)
mod <- fit(mod)
ans_obtained <- make_linpred_raw(mod, point = TRUE)
m1 <- make_combined_matrix_effect_outcome(mod)
m2 <- make_combined_matrix_effectfree_effect(mod)
ans_expected <- as.double(m1 %*% m2 %*% mod$point_effectfree)
expect_equal(ans_obtained, ans_expected)
})
## 'make_point_est_effects' ---------------------------------------------------
test_that("'make_point_est_effects' works with valid inputs", {
set.seed(0)
data <- expand.grid(age = 0:9,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(age ~ RW(sd = 0)) |>
set_prior(age:sex ~ RW(sd = 0))
mod <- fit(mod)
ans_obtained <- make_point_est_effects(mod)
int <- unname(mod$point_effectfree[1])
age <- c(0, unname(mod$point_effectfree[2:10]))
sex <- unname(mod$point_effectfree[11:12])
agesex <- c(0, unname(mod$point_effectfree[13:21]),
0, unname(mod$point_effectfree[22:30]))
ans_expected <- list("(Intercept)" = int,
age = age,
sex = sex,
"age:sex" = agesex)
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_point_est_effects' throws correct error when not fitted", {
set.seed(0)
data <- expand.grid(age = 0:9,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
expect_error(make_point_est_effects(mod),
"Internal error: Model not fitted.")
})
## 'make_scaled_eigen' --------------------------------------------------------
## See also tests for rvnorm_eigen
test_that("'make_scaled_eigen' works with positive definite matrix", {
set.seed(0)
prec <- solve(crossprod(matrix(rnorm(25), 5)))
ans <- make_scaled_eigen(prec)
expect_identical(dim(ans), dim(prec))
})
test_that("'make_scaled_eigen' works with non-negative definite matrix", {
set.seed(0)
prec <- solve(crossprod(matrix(rnorm(25), 5)))
prec[5,] <- 0
prec[,5] <- 0
ans <- make_scaled_eigen(prec)
expect_identical(dim(ans), dim(prec))
})
## 'make_spline' --------------------------------------------------------------
test_that("'make_spline' works", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age ~ RW(sd = 0))
mod <- set_prior(mod, age:sex ~ Sp(n = 5))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_spline(mod = mod, effectfree = effectfree)
ans_expected <- effectfree[23:32,]
expect_equal(ans_obtained, ans_expected)
})
## 'make_svd' -----------------------------------------------------------------
test_that("'make_svd' works - SVD_RW, random", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:sex ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD, n_comp = 5))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_svd(mod = mod, effectfree = effectfree)
ans_expected <- effectfree[24:63,]
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_svd' works - SVD_RW, zero", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:sex ~ SVD(HMD, n_comp = 5))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD, n_comp = 5, sd = 0))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_svd(mod = mod, effectfree = effectfree)
ans_expected <- rbind(effectfree[24:33,],
matrix(0, nrow = 5, ncol = 5),
effectfree[34:58,])
expect_equal(ans_obtained, ans_expected)
})
test_that("'make_svd' works - SVD_AR", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 20)
formula <- deaths ~ age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, age:time ~ SVD_AR1(HMD, n_comp = 1))
mod <- set_n_draw(mod, n = 5)
mod <- fit(mod)
effectfree <- mod$draws_effectfree
ans_obtained <- make_svd(mod = mod, effectfree = effectfree)
ans_expected <- effectfree[22:27,]
expect_equal(ans_obtained, ans_expected)
})
## 'make_term_components' -----------------------------------------------------
test_that("'make_term_components' works - no disp", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_disp(mod, mean = 0)
mod <- set_prior(mod, age ~ Sp())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_term_components(mod)
expect_identical(length(ans), length(comp))
})
test_that("'make_term_components' works - has hyperrand", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_term_components(mod)
expect_identical(length(ans), length(comp))
})
test_that("'make_term_components' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 1)
mod <- fit(mod)
comp <- make_comp_components(mod)
ans <- make_term_components(mod)
expect_identical(length(ans), length(comp))
})
## 'make_term_spline' ------------------------------------------------------------
test_that("'make_term_spline' works - no spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
ans_obtained <- make_term_spline(mod)
ans_expected <- factor()
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_term_spline' works - has spline", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ Sp(n_comp = 5))
ans_obtained <- make_term_spline(mod)
ans_expected <- factor(rep("age", times = 5))
expect_identical(ans_obtained, ans_expected)
})
## 'make_term_svd' ------------------------------------------------------------
test_that("'make_term_svd' works - no svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 3)
ans_obtained <- make_term_svd(mod)
ans_expected <- factor()
expect_identical(ans_obtained, ans_expected)
})
test_that("'make_term_svd' works - has svd", {
set.seed(0)
data <- expand.grid(age = poputils::age_labels(type = "lt", max = 60),
time = 2000:2005,
sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age * sex + age * time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_n_draw(mod, n = 5)
mod <- set_prior(mod, age ~ SVD(HMD))
mod <- set_prior(mod, age:time ~ SVD_RW(HMD))
ans_obtained <- make_term_svd(mod)
ans_expected <- factor(c(rep("age", times = 3),
rep("age:time", times = 3 * 6)))
expect_identical(ans_obtained, ans_expected)
})
## 'make_terms_hyperrand' ----------------------------------------------------
test_that("'make_terms_hyperrand' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex:time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- set_prior(mod, sex:time ~ Lin())
ans_obtained <- make_terms_hyperrand(mod)
ans_expected <- factor(rep("sex:time", 26),
levels = c("(Intercept)", "age", "sex:time"))
expect_identical(ans_obtained, ans_expected)
})
## 'make_unconstr_dimnames_by' ------------------------------------------------
test_that("'make_unconstr_dimnames_by' works", {
dimnames_term <- list(age = 1:5, time = 1:6, reg = 1:3)
ans_obtained <- make_unconstr_dimnames_by(i_along = 2L,
dimnames_term = dimnames_term)
ans_expected <- list(age = paste0("age", 1:4), reg = paste0("reg", 1:2))
expect_identical(ans_obtained, ans_expected)
})
## 'paste_dot' ----------------------------------------------------------------
test_that("'paste_dot' works with valid inputs", {
expect_identical(paste_dot(1:3, 3:1), c("1.3", "2.2", "3.1"))
})
## 'rescale_lin_intercept' ----------------------------------------------------
test_that("'rescale_lin_intercept' works with valid inputs", {
set.seed(0)
slope <- rvec::rnorm_rvec(n = 4, n_draw = 10)
effect <- rvec::rnorm_rvec(n = 20, n_draw = 10)
matrix_along_by <- matrix(0:19, nrow = 5)
ans_obtained <- rescale_lin_intercept(slope = slope,
effect = effect,
matrix_along_by = matrix_along_by)
ans_expected <- slope
for (i in 1:4)
ans_expected[i] <- mean(effect[1:5 + (i - 1) * 5]) - mean(slope[i] * (1:5))
expect_equal(ans_obtained, ans_expected)
})
## 'make_transforms_hyper' ----------------------------------------------------
test_that("'make_transforms_hyper' works", {
set.seed(0)
data <- expand.grid(age = 0:9, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ age + sex + time
mod <- mod_pois(formula = formula,
data = data,
exposure = popn)
mod <- fit(mod)
ans_obtained <- make_transforms_hyper(mod)
ans_expected <- list(exp, exp)
expect_identical(unname(ans_obtained), ans_expected,
ignore_function_env = TRUE)
})
## 'infer_trend_cyc_seas_err_forecast' --------------------------------------------
test_that("'infer_trend_cyc_seas_err_forecast' works", {
set.seed(0)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * time + age
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(time ~ RW2()) |>
set_prior(sex:time ~ Lin()) |>
fit()
mod <- set_n_draw(mod, 5)
mod <- fit(mod)
comp_est <- components(mod)
comp_forecast <- forecast(mod, labels = 2006:2007, output = "components")
dimnames_terms_forecast <- make_dimnames_terms_forecast(dimnames_terms = mod$dimnames_terms,
var_time = mod$var_time,
labels_forecast = 2006:2007,
time_only = TRUE)
ans_obtained <- infer_trend_cyc_seas_err_forecast(components = comp_forecast,
priors = mod$priors,
dimnames_terms = dimnames_terms_forecast,
var_time = mod$var_time,
var_age = mod$var_age)
expect_equal(ans_obtained[1:3], comp_forecast[1:3])
})
## 'infer_trend_cyc_seas_err_seasfix_forecast' --------------------------------
test_that("'infer_trend_cyc_seas_err_seasfix_forecast' works", {
set.seed(0)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * time + age
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ RW_Seas(n = 3)) |>
set_n_draw(n = 10) |>
fit()
components <- components(mod)
ans <- infer_trend_cyc_seas_err_seasfix_forecast(prior = mod$priors[["sex:time"]],
dimnames_term = mod$dimnames_terms[["sex:time"]],
var_time = mod$var_time,
var_age = mod$var_age,
components = components)
season <- ans$.fitted[ans$term == "sex:time" & ans$component == "season"]
trend <- ans$.fitted[ans$term == "sex:time" & ans$component == "trend"]
effect <- ans$.fitted[ans$term == "sex:time" & ans$component == "effect"]
expect_equal(effect, season + trend)
})
## 'infer_trend_cyc_seas_err_seasvary_forecast' -------------------------------
test_that("'infer_trend_cyc_seas_err_seasvary_forecast' works", {
set.seed(0)
data <- expand.grid(age = 0:4, time = 2000:2005, sex = c("F", "M"))
data$popn <- rpois(n = nrow(data), lambda = 100)
data$deaths <- rpois(n = nrow(data), lambda = 10)
formula <- deaths ~ sex * time + age
mod <- mod_pois(formula = formula,
data = data,
exposure = popn) |>
set_prior(sex:time ~ RW_Seas(n = 3, s = 1)) |>
set_n_draw(n = 10) |>
fit()
components <- components(mod)
ans <- infer_trend_cyc_seas_err_seasvary_forecast(prior = mod$priors[["sex:time"]],
dimnames_term = mod$dimnames_terms[["sex:time"]],
var_time = mod$var_time,
var_age = mod$var_age,
components = components)
season <- ans$.fitted[ans$term == "sex:time" & ans$component == "season"]
trend <- ans$.fitted[ans$term == "sex:time" & ans$component == "trend"]
effect <- ans$.fitted[ans$term == "sex:time" & ans$component == "effect"]
expect_equal(effect, season + trend)
})
## 'rvec_to_mean' -------------------------------------------------------------
test_that("'rvec_to_mean' works with valid inputs", {
data <- tibble::tibble(a = 1:3,
b = rvec::rvec(matrix(1:12, nr = 3)),
c = rvec::rvec(matrix(FALSE, nr = 3, nc = 2)),
d = c("a", "b", "c"))
ans_obtained <- rvec_to_mean(data)
ans_expected <- tibble::tibble(a = 1:3,
b = rowMeans(matrix(1:12, nr = 3)),
c = rowMeans(matrix(FALSE, nr = 3, nc = 2)),
d = c("a", "b", "c"))
expect_identical(ans_obtained, ans_expected)
})
## 'rmvnorm_chol', 'rmvnorm_eigen' --------------------------------------------
test_that("'rmvnorm_chol' and 'rmvnorm_eigen' give the same answer", {
set.seed(0)
prec <- crossprod(matrix(rnorm(25), 5))
mean <- rnorm(5)
R_prec <- chol(prec)
scaled_eigen <- make_scaled_eigen(prec)
ans_chol <- rmvnorm_chol(n = 100000, mean = mean, R_prec = R_prec)
ans_eigen <- rmvnorm_eigen(n = 100000, mean = mean, scaled_eigen = scaled_eigen)
expect_equal(rowMeans(ans_chol), rowMeans(ans_eigen), tolerance = 0.02)
expect_equal(cov(t(ans_chol)), cov(t(ans_eigen)), tolerance = 0.02)
})
## 'sort_components' ----------------------------------------------------------
test_that("'sort_components' works with valid inputs", {
components <- tibble::tribble(~term, ~component, ~level,
"sex", "effect", "m",
"time", "season", "2000",
"sex", "hyper", "sd",
"sex", "effect", "f",
"(Intercept)", "effect", "(Intercept)",
"time", "effect", "2000")
mod <- list(formula = deaths ~ time + sex)
ans_obtained <- sort_components(components = components, mod = mod)
ans_expected <- tibble::tribble(~term, ~component, ~level,
"(Intercept)", "effect", "(Intercept)",
"time", "effect", "2000",
"time", "season", "2000",
"sex", "effect", "m",
"sex", "effect", "f",
"sex", "hyper", "sd")
expect_identical(ans_obtained, ans_expected)
})
test_that("'sort_components' raises correct effor with invalid component", {
components <- tibble::tribble(~term, ~component, ~level,
"(Intercept)", "effect", "(Intercept)",
"time", "season", "2000",
"sex", "wrong", "sd")
expect_error(sort_components(components, mod = list(formula = deaths ~ age)),
"Internal error: \"wrong\" not a valid value for `component`.")
})
## 'transform_hyper_ar' -------------------------------------------------------
test_that("'transform_hyper_ar' works with 'bage_prior_ar - AR1'", {
shifted_invlogit <- function(x) {
ans <- exp(x) / (1 + exp(x))
0.18 * ans + 0.8
}
l <- transform_hyper_ar(prior = AR1())
expect_equal(l[[1]](0.35), shifted_invlogit(0.35))
expect_equal(l[[2]](0.35), exp(0.35))
})
test_that("'transform_hyper_ar' works with 'bage_prior_svd_ar - AR'", {
shifted_invlogit <- function(x) {
ans <- exp(x) / (1 + exp(x))
2 * ans - 1
}
l <- transform_hyper_ar(prior = SVD_AR(HMD))
expect_equal(l[[1]](0.35), shifted_invlogit(0.35))
expect_equal(l[[2]](0.35), shifted_invlogit(0.35))
expect_equal(l[[3]](0.35), exp(0.35))
})
Try the bage package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.