context("project()")
# object, predictor_terms, ndraws / nclusters, nterms ---------------------
test_that(paste(
"`object` of class `refmodel` and arguments `predictor_terms` and `ndraws`",
"(or `nclusters`) work"
), {
skip_if_not(run_prj)
for (tstsetup in names(prjs)) {
args_prj_i <- args_prj[[tstsetup]]
ndr_ncl <- ndr_ncl_dtls(args_prj_i)
projection_tester(prjs[[tstsetup]],
refmod_expected = refmods[[args_prj_i$tstsetup_ref]],
prd_trms_expected = args_prj_i$predictor_terms,
nprjdraws_expected = ndr_ncl$nprjdraws,
with_clusters = ndr_ncl$clust_used,
info_str = tstsetup)
}
})
test_that("invalid `predictor_terms` warns or fails", {
skip_if_not(run_prj)
tstsetups <- grep("\\.glm\\.gauss.*\\.prd_trms_x\\.clust$", names(prjs),
value = TRUE)
for (tstsetup in tstsetups) {
args_prj_i <- args_prj[[tstsetup]]
# Non-`vsel` object combined with `predictor_terms = NULL`:
expect_error(
do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]],
predictor_terms = NULL),
excl_nonargs(args_prj_i, nms_excl_add = "predictor_terms")
)),
paste("^Please provide an `object` of class `vsel` or use argument",
"`predictor_terms`\\.$"),
info = tstsetup
)
# Repeating `predictor_terms`:
p_long <- do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]],
predictor_terms = rep_len(args_prj_i$predictor_terms,
length.out = 1e4)),
excl_nonargs(args_prj_i, nms_excl_add = "predictor_terms")
))
expect_equal(p_long, prjs[[tstsetup]], info = tstsetup)
# Invalid type:
for (prd_trms_crr in list(2, 1:3, list(prd_trms_x, prd_trms_x))) {
tstsetup_crr <- paste(tstsetup, paste(prd_trms_crr, collapse = ","),
sep = "__")
expect_error(
do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]],
predictor_terms = prd_trms_crr),
excl_nonargs(args_prj_i, nms_excl_add = "predictor_terms")
)),
paste(
"^is\\.null\\(predictor_terms\\) \\|\\|",
"is\\.vector\\(predictor_terms, \"character\"\\) is not TRUE$"
),
info = tstsetup_crr
)
}
# Should be working, but result in a projection onto the intercept-only
# submodel:
for (prd_trms_crr in list("1",
c("some_dummy_string", "another_dummy_string"))) {
tstsetup_crr <- paste(tstsetup, paste(prd_trms_crr, collapse = ","),
sep = "__")
expect_warning(
p <- do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]],
predictor_terms = prd_trms_crr),
excl_nonargs(args_prj_i, nms_excl_add = "predictor_terms")
)),
paste("The following element\\(s\\) of `predictor_terms` could not be",
"found in the table of possible predictor terms"),
info = tstsetup_crr
)
projection_tester(p,
refmod_expected = refmods[[args_prj_i$tstsetup_ref]],
prd_trms_expected = character(),
nprjdraws_expected = nclusters_pred_tst,
with_clusters = TRUE,
info_str = tstsetup_crr)
}
}
})
test_that("`object` of class `stanreg` or `brmsfit` works", {
skip_if_not(run_prj)
tstsetups <- grep("\\.brnll\\..*\\.prd_trms_x\\.clust$", names(prjs),
value = TRUE)
for (tstsetup in tstsetups) {
args_prj_i <- args_prj[[tstsetup]]
p_fit <- do.call(project, c(
list(object = fits[[args_prj_i$tstsetup_fit]]),
excl_nonargs(args_prj_i),
excl_nonargs(args_ref[[args_prj_i$tstsetup_ref]])
))
expect_identical(p_fit, prjs[[tstsetup]], ignore.environment = TRUE,
info = tstsetup)
}
})
test_that(paste(
"`object` of class `vsel` (created by varsel()) and arguments `nclusters`",
"and `nterms` work"
), {
skip_if_not(run_vs)
for (tstsetup in names(prjs_vs)) {
tstsetup_vs <- args_prj_vs[[tstsetup]]$tstsetup_vsel
stopifnot(length(tstsetup_vs) > 0)
mod_crr <- args_vs[[tstsetup_vs]]$mod_nm
fam_crr <- args_vs[[tstsetup_vs]]$fam_nm
nterms_crr <- args_prj_vs[[tstsetup]]$nterms
if (is.null(nterms_crr)) {
nterms_crr <- suggest_size(vss[[tstsetup_vs]], warnings = FALSE)
}
if (length(nterms_crr) == 1) {
prd_trms_expected_crr <- vss[[tstsetup_vs]]$predictor_ranking[
seq_len(nterms_crr)
]
projection_tester(
prjs_vs[[tstsetup]],
refmod_expected = refmods[[args_prj_vs[[tstsetup]]$tstsetup_ref]],
prd_trms_expected = prd_trms_expected_crr,
nprjdraws_expected = args_prj_vs[[tstsetup]]$nclusters,
with_clusters = TRUE,
info_str = tstsetup
)
# Check that projecting from the `vsel` object onto a single submodel
# gives the same output as projecting the reference model onto that
# submodel directly:
tstsetup_tries <- grep(
paste0("^",
gsub("\\.", "\\\\.",
sub("(with_offs|without_offs).*", "\\1", tstsetup)),
".*\\.clust$"),
names(prjs),
value = TRUE
)
match_prj <- sapply(tstsetup_tries, function(tstsetup_try) {
setequal(prd_trms_expected_crr, prjs[[tstsetup_try]]$predictor_terms) &&
args_prj_vs[[tstsetup]]$prj_nm == args_prj[[tstsetup_try]]$prj_nm
})
tstsetup_match_prj <- tstsetup_tries[match_prj]
if (length(tstsetup_match_prj) == 1) {
if (identical(prjs_vs[[tstsetup]]$predictor_terms,
prjs[[tstsetup_match_prj]]$predictor_terms)) {
expect_identical(prjs_vs[[tstsetup]], prjs[[tstsetup_match_prj]],
ignore.environment = TRUE, info = tstsetup)
} else {
expect_setequal(prjs_vs[[tstsetup]]$predictor_terms,
prjs[[tstsetup_match_prj]]$predictor_terms)
expect_equal(
lapply(seq_along(prjs_vs[[tstsetup]]$outdmin), function(s_idx) {
outdmin_s <- prjs_vs[[tstsetup]]$outdmin[[s_idx]]
outdmin_s$beta <- outdmin_s$beta[
rownames(prjs[[tstsetup_match_prj]]$outdmin[[s_idx]]$beta), ,
drop = FALSE
]
outdmin_s$formula <- prjs[[tstsetup_match_prj]]$outdmin[[s_idx]][[
"formula"
]]
outdmin_s$x <- outdmin_s$x[
, colnames(prjs[[tstsetup_match_prj]]$outdmin[[s_idx]]$x),
drop = FALSE
]
return(outdmin_s)
}),
prjs[[tstsetup_match_prj]]$outdmin,
tolerance = 1e1 * .Machine$double.eps, info = tstsetup
)
prj_nms <- names(prjs_vs[[tstsetup]])
expect_identical(prj_nms, names(prjs[[tstsetup_match_prj]]),
info = tstsetup)
prj_el_excl <- !prj_nms %in% c("predictor_terms", "outdmin")
expect_equal(prjs_vs[[tstsetup]][prj_el_excl],
prjs[[tstsetup_match_prj]][prj_el_excl],
tolerance = .Machine$double.eps, info = tstsetup)
}
} else if (length(tstsetup_match_prj) > 1) {
stop("Unexpected number of matches.")
}
} else {
proj_list_tester(
prjs_vs[[tstsetup]],
len_expected = length(nterms_crr),
is_seq = all(diff(nterms_crr) == 1),
info_str = tstsetup,
refmod_expected = refmods[[args_prj_vs[[tstsetup]]$tstsetup_ref]],
nprjdraws_expected = args_prj_vs[[tstsetup]]$nclusters,
with_clusters = TRUE,
prjdraw_weights_expected = prjs_vs[[tstsetup]][[1]]$wdraws_prj
)
}
}
})
test_that(paste(
"`object` of class `vsel` (created by cv_varsel()) and arguments",
"`nclusters` and `nterms` work"
), {
skip_if_not(run_cvvs)
for (tstsetup in names(prjs_cvvs)) {
tstsetup_cvvs <- args_prj_cvvs[[tstsetup]]$tstsetup_vsel
stopifnot(length(tstsetup_cvvs) > 0)
mod_crr <- args_cvvs[[tstsetup_cvvs]]$mod_nm
fam_crr <- args_cvvs[[tstsetup_cvvs]]$fam_nm
nterms_crr <- args_prj_cvvs[[tstsetup]]$nterms
if (is.null(nterms_crr)) {
nterms_crr <- suggest_size(cvvss[[tstsetup_cvvs]], warnings = FALSE)
}
if (length(nterms_crr) == 1) {
prd_trms_expected_crr <- cvvss[[tstsetup_cvvs]]$predictor_ranking[
seq_len(nterms_crr)
]
projection_tester(
prjs_cvvs[[tstsetup]],
refmod_expected = refmods[[args_prj_cvvs[[tstsetup]]$tstsetup_ref]],
prd_trms_expected = prd_trms_expected_crr,
nprjdraws_expected = args_prj_cvvs[[tstsetup]]$nclusters,
with_clusters = TRUE,
info_str = tstsetup
)
# Check that projecting from the `vsel` object onto a single submodel
# gives the same output as projecting the reference model onto that
# submodel directly:
tstsetup_tries <- grep(
paste0("^",
gsub("\\.", "\\\\.",
sub("(with_offs|without_offs).*", "\\1", tstsetup)),
".*\\.clust$"),
names(prjs),
value = TRUE
)
match_prj <- sapply(tstsetup_tries, function(tstsetup_try) {
setequal(prd_trms_expected_crr, prjs[[tstsetup_try]]$predictor_terms) &&
args_prj_cvvs[[tstsetup]]$prj_nm == args_prj[[tstsetup_try]]$prj_nm
})
tstsetup_match_prj <- tstsetup_tries[match_prj]
if (length(tstsetup_match_prj) == 1) {
if (identical(prjs_cvvs[[tstsetup]]$predictor_terms,
prjs[[tstsetup_match_prj]]$predictor_terms)) {
expect_identical(prjs_cvvs[[tstsetup]], prjs[[tstsetup_match_prj]],
ignore.environment = TRUE, info = tstsetup)
} else {
expect_setequal(prjs_cvvs[[tstsetup]]$predictor_terms,
prjs[[tstsetup_match_prj]]$predictor_terms)
expect_equal(
lapply(seq_along(prjs_cvvs[[tstsetup]]$outdmin), function(s_idx) {
outdmin_s <- prjs_cvvs[[tstsetup]]$outdmin[[s_idx]]
outdmin_s$beta <- outdmin_s$beta[
rownames(prjs[[tstsetup_match_prj]]$outdmin[[s_idx]]$beta), ,
drop = FALSE
]
outdmin_s$formula <- prjs[[tstsetup_match_prj]]$outdmin[[s_idx]][[
"formula"
]]
outdmin_s$x <- outdmin_s$x[
, colnames(prjs[[tstsetup_match_prj]]$outdmin[[s_idx]]$x),
drop = FALSE
]
return(outdmin_s)
}),
prjs[[tstsetup_match_prj]]$outdmin,
tolerance = 1e1 * .Machine$double.eps, info = tstsetup
)
prj_nms <- names(prjs_cvvs[[tstsetup]])
expect_identical(prj_nms, names(prjs[[tstsetup_match_prj]]),
info = tstsetup)
prj_el_excl <- !prj_nms %in% c("predictor_terms", "outdmin")
expect_equal(prjs_cvvs[[tstsetup]][prj_el_excl],
prjs[[tstsetup_match_prj]][prj_el_excl],
tolerance = .Machine$double.eps, info = tstsetup)
}
} else if (length(tstsetup_match_prj) > 1) {
stop("Unexpected number of matches.")
}
} else {
proj_list_tester(
prjs_cvvs[[tstsetup]],
len_expected = length(nterms_crr),
is_seq = all(diff(nterms_crr) == 1),
info_str = tstsetup,
refmod_expected = refmods[[args_prj_cvvs[[tstsetup]]$tstsetup_ref]],
nprjdraws_expected = args_prj_cvvs[[tstsetup]]$nclusters,
with_clusters = TRUE,
prjdraw_weights_expected = prjs_cvvs[[tstsetup]][[1]]$wdraws_prj
)
}
}
})
test_that("`object` not of class `vsel` and missing `predictor_terms` fails", {
skip_if_not(length(fits) > 0)
expect_error(
project(fits[[1]]),
paste("^Please provide an `object` of class `vsel` or use argument",
"`predictor_terms`\\.$")
)
expect_error(
project(refmods[[1]]),
paste("^Please provide an `object` of class `vsel` or use argument",
"`predictor_terms`\\.$")
)
})
# nterms ------------------------------------------------------------------
test_that("invalid `nterms` fails", {
skip_if_not(run_vs)
tstsetups <- head(grep("\\.glm\\.gauss", names(vss), value = TRUE), 1)
for (tstsetup in tstsetups) {
for (nterms_crr in nterms_unavail) {
expect_error(project(vss[[!!tstsetup]], nterms = !!nterms_crr),
paste("Cannot perform the projection with", max(nterms_crr),
"variables"))
}
for (nterms_crr in list(neg = -1, char = "a", dafr = dat)) {
expect_error(project(vss[[!!tstsetup]], nterms = !!nterms_crr),
"must contain non-negative values")
}
}
})
# seed --------------------------------------------------------------------
test_that("non-clustered projection does not require a seed", {
skip_if_not(run_prj)
# This test is important to ensure that we don't have to set a seed where we
# don't expect it to be necessary.
tstsetups <- grep("\\.noclust$|\\.default_ndr_ncl$", names(prjs),
value = TRUE)
for (tstsetup in tstsetups) {
args_prj_i <- args_prj[[tstsetup]]
p_orig <- prjs[[tstsetup]]
rand_new1 <- runif(1) # Just to advance `.Random.seed[2]`.
# Use suppressWarnings() because test_that() somehow redirects stderr() and
# so throws warnings that projpred wants to capture internally:
p_new <- suppressWarnings(do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]]),
excl_nonargs(args_prj_i, nms_excl_add = "seed")
)))
if (args_prj_i$mod_nm %in% c("glmm", "gamm") &&
any(grepl("\\|", args_prj_i$predictor_terms))) {
if (getOption("projpred.mlvl_pred_new", FALSE)) {
# In this case, the multilevel submodel fitters (fit_glmer_callback(),
# fit_gamm_callback(), fit_cumul_mlvl(), fit_categ_mlvl()) should still
# be deterministic, but the prediction from the fitted submodels is not
# (because of the group-level effects drawn randomly by repair_re() (for
# all group levels; here, only the existing ones are relevant)). Thus,
# we cannot test the whole project() output, but need to restrict
# ourselves to the output of as.matrix.projection().
if (args_prj_i$mod_nm == "gamm") {
# Skipping GAMMs because of issue #131.
# TODO (GAMMs): Fix this.
next
}
prjmat_orig <- as.matrix(p_orig)
prjmat_new <- as.matrix(p_new)
if (args_prj_i$fam_nm == "gauss" || args_prj_i$prj_nm == "latent") {
# The projected dispersion parameter is affected by the group-level
# effects drawn randomly by repair_re() (for all group levels):
prjmat_new[, "sigma"] <- prjmat_orig[, "sigma"]
}
expect_equal(prjmat_new, prjmat_orig, info = tstsetup,
tolerance = .Machine$double.eps)
# To facilitate the `if` conditions here:
p_new <- p_orig
} else if (args_prj_i$prj_nm == "augdat" &&
args_prj_i$fam_nm == "cumul" && args_prj_i$mod_nm == "glmm") {
for (idx_s in seq_along(p_new$outdmin)) {
if (!is.null(p_new$outdmin[[idx_s]][["L"]])) {
# We could also use `sparseMatrix` instead of `Matrix`:
expect_equal(as(p_new$outdmin[[idx_s]][["L"]], "Matrix"),
as(p_orig$outdmin[[idx_s]][["L"]], "Matrix"),
info = tstsetup)
p_new$outdmin[[idx_s]][["L"]] <- p_orig$outdmin[[idx_s]][["L"]]
}
}
}
}
expect_equal(p_new, p_orig, info = tstsetup)
}
})
test_that("`seed` works (and restores the RNG state afterwards)", {
skip_if_not(run_prj)
tstsetups <- grep("\\.glm\\.gauss.*\\.prd_trms_x\\.clust$", names(prjs),
value = TRUE)
for (tstsetup in tstsetups) {
args_prj_i <- args_prj[[tstsetup]]
p_orig <- prjs[[tstsetup]]
rand_orig <- runif(1) # Just to advance `.Random.seed[2]`.
.Random.seed_new1 <- .Random.seed
p_new <- do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]],
seed = args_prj_i$seed + 10L),
excl_nonargs(args_prj_i, nms_excl_add = "seed")
))
.Random.seed_new2 <- .Random.seed
rand_new <- runif(1) # Just to advance `.Random.seed[2]`.
.Random.seed_repr1 <- .Random.seed
p_repr <- do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]]),
excl_nonargs(args_prj_i)
))
.Random.seed_repr2 <- .Random.seed
# Expected equality:
expect_equal(p_repr, p_orig, info = tstsetup)
expect_equal(.Random.seed_new2, .Random.seed_new1, info = tstsetup)
expect_equal(.Random.seed_repr2, .Random.seed_repr1, info = tstsetup)
# Expected inequality:
expect_false(isTRUE(all.equal(p_new, p_orig)), info = tstsetup)
expect_false(isTRUE(all.equal(rand_new, rand_orig)), info = tstsetup)
expect_false(isTRUE(all.equal(.Random.seed_repr2, .Random.seed_new2)),
info = tstsetup)
}
})
# regul -------------------------------------------------------------------
test_that("for GLMs, `regul` has an expected effect", {
skip_if_not(run_prj)
regul_tst <- c(regul_default, 1e-1, 1e2)
stopifnot(regul_tst[1] == regul_default)
stopifnot(all(diff(regul_tst) > 0))
tstsetups <- grep("\\.glm\\..*\\.clust$", names(prjs), value = TRUE)
tstsetups <- grep(fam_nms_aug_regex, tstsetups, value = TRUE, invert = TRUE)
for (tstsetup in tstsetups) {
args_prj_i <- args_prj[[tstsetup]]
ndr_ncl <- ndr_ncl_dtls(args_prj_i)
# Calculate the objects for which to run checks:
ssq_regul_alpha <- rep(NA, length(regul_tst))
ssq_regul_beta <- rep(NA, length(regul_tst))
for (j in seq_along(regul_tst)) {
# Run project() if necessary:
if (regul_tst[j] == regul_default) {
prj_regul <- prjs[[tstsetup]]
} else {
prj_regul <- do.call(project, c(
list(object = refmods[[args_prj_i$tstsetup_ref]],
regul = regul_tst[j]),
excl_nonargs(args_prj_i)
))
projection_tester(
prj_regul,
refmod_expected = refmods[[args_prj_i$tstsetup_ref]],
prd_trms_expected = args_prj_i$predictor_terms,
nprjdraws_expected = ndr_ncl$nprjdraws,
with_clusters = ndr_ncl$clust_used,
info_str = paste(tstsetup, j, sep = "__")
)
}
# Run as.matrix.projection():
prjmat <- as.matrix(prj_regul, nm_scheme = "brms",
allow_nonconst_wdraws_prj = ndr_ncl$clust_used_gt1)
# Reduce to only those columns which are necessary here:
prjmat <- prjmat[, grep("^b_", colnames(prjmat)), drop = FALSE]
# Posterior means:
prjmat_mean <- colMeans(prjmat)
# Calculate the Euclidean norm for intercept and coefficients:
ssq_regul_alpha[j] <- prjmat_mean["b_Intercept"]^2
coef_colnms <- setdiff(names(prjmat_mean), "b_Intercept")
if (length(coef_colnms) > 0) {
ssq_regul_beta[j] <- sum(prjmat_mean[coef_colnms]^2)
}
}
# Checks:
if (length(args_prj_i$predictor_terms) == 0) {
# For an intercept-only model:
expect_length(unique(ssq_regul_alpha), 1)
stopifnot(all(is.na(ssq_regul_beta)))
} else {
# All other (i.e., not intercept-only) models (note: as discussed at issue
# #169, the intercept is not tested here to stay the same):
for (j in seq_along(ssq_regul_beta)[-1]) {
expect_lt(ssq_regul_beta[!!j], ssq_regul_beta[j - 1])
}
}
}
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.