Skip to content

Commit

Permalink
can't yet handle id in smooth terms
Browse files Browse the repository at this point in the history
  • Loading branch information
Nicholas Clark committed May 16, 2024
1 parent b434414 commit 9ced215
Show file tree
Hide file tree
Showing 4 changed files with 135 additions and 117 deletions.
11 changes: 10 additions & 1 deletion R/mvgam.R
Original file line number Diff line number Diff line change
Expand Up @@ -854,9 +854,18 @@ mvgam = function(formula,
# For monotonic smooths, need to determine which direction to place
# coefficient constraints
smooth_labs <- do.call(rbind, lapply(seq_along(ss_gam$smooth), function(x){
data.frame(label = ss_gam$smooth[[x]]$label, class = class(ss_gam$smooth[[x]])[1])
data.frame(label = ss_gam$smooth[[x]]$label,
class = class(ss_gam$smooth[[x]])[1],
id = ifelse(is.null(ss_gam$smooth[[x]]$id),
NA, ss_gam$smooth[[x]]$id))
}))

# Check for 'id' arguments, which are not yet supported
if(any(!is.na(smooth_labs$id))){
stop('smooth terms with the "id" argument not yet supported by mvgam',
call. = FALSE)
}

if(any(smooth_labs$class == 'random.effect')){
re_smooths <- smooth_labs %>%
dplyr::filter(class == 'random.effect') %>%
Expand Down
Binary file modified src/mvgam.dll
Binary file not shown.
Binary file modified tests/testthat/Rplots.pdf
Binary file not shown.
241 changes: 125 additions & 116 deletions tests/testthat/test-mvgam.R
Original file line number Diff line number Diff line change
@@ -1,5 +1,130 @@
context("mvgam")

test_that("family must be correctly specified", {
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = beta_data$data_train,
family = 'besta',
run_model = FALSE),
'family not recognized')
})

test_that("response variable must be specified", {
expect_error(mod <- mvgam( ~ s(season),
trend_model = 'AR1',
data = beta_data$data_train,
family = betar(),
run_model = FALSE),
'response variable is missing from formula')
})

test_that("id to link smooths not allowed yet", {
expect_error(mod <- mvgam(y ~ s(time, id = 1) +
s(time, by = series, id = 1),
data = beta_data$data_train,
family = betar(),
run_model = FALSE),
'smooth terms with the "id" argument not yet supported by mvgam')
})

test_that("response variable must follow family-specific restrictions", {
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = gaus_data$data_train,
family = lognormal(),
run_model = FALSE),
'Values <= 0 not allowed for lognormal responses')

expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = gaus_data$data_train,
family = poisson(),
run_model = FALSE),
'Values < 0 not allowed for count family responses')
})

test_that("trend_model must be correctly specified", {
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR11',
data = beta_data$data_train,
family = betar(),
run_model = FALSE))
})

test_that("outcome variable must be present in data", {
data = data.frame(out = rnorm(100),
temp = rnorm(100),
time = 1:100)
expect_error(mod <- mvgam(formula = y ~ dynamic(temp, rho = 20),
data = data,
family = gaussian(),
run_model = FALSE),
'variable y not found in data')
})

test_that("series levels must match unique entries in series", {
levels(beta_data$data_train$series) <- paste0('series_', 1:6)
expect_error(mvgam(y ~ s(season),
trend_model = 'GP',
data = beta_data$data_train,
newdata = beta_data$data_test,
family = betar(),
run_model = FALSE))
})

test_that("missing values not allowed in predictors", {
# Include missing vals in training data
simdat <- sim_mvgam()
simdat$data_train$season[4] <- NA
expect_error(mvgam(y ~ s(season),
trend_model = 'GP',
data = simdat$data_train,
newdata = simdat$data_test,
run_model = FALSE))

# Include missing vals in testing data
simdat <- sim_mvgam()
simdat$data_test$season[4] <- NA
expect_error(mvgam(y ~ s(season),
data = simdat$data_train,
newdata = simdat$data_test,
run_model = FALSE))
})

test_that("all series must have observations for all unique timepoints", {
data <- sim_mvgam()
data$data_train <- data$data_train[-2,]
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = data$data_train,
family = poisson(),
run_model = FALSE),
'One or more series in data is missing observations for one or more timepoints')

data <- sim_mvgam()
data$data_test <- data$data_test[-2,]
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = data$data_train,
newdata = data$data_test,
family = poisson(),
run_model = FALSE),
'One or more series in newdata is missing observations for one or more timepoints')
})

test_that("rho argument must be positive numeric", {
data = data.frame(out = rnorm(100),
temp = rnorm(100),
time = 1:100)
expect_error(mod <- mvgam(formula = out ~ dynamic(temp, rho = -1),
data = data,
family = gaussian(),
run_model = FALSE),
'Argument "rho" in dynamic() must be a positive value',
fixed = TRUE)
})

# Skip remaining tests on CRAN as they are slightly time-consuming
skip_on_cran()

test_that("JAGS setups should work", {
Expand Down Expand Up @@ -82,24 +207,6 @@ test_that("JAGS setups should work", {
expect_true(inherits(mod, 'mvgam_prefit'))
})

test_that("family must be correctly specified", {
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = beta_data$data_train,
family = 'besta',
run_model = FALSE),
'family not recognized')
})

test_that("response variable must be specified", {
expect_error(mod <- mvgam( ~ s(season),
trend_model = 'AR1',
data = beta_data$data_train,
family = betar(),
run_model = FALSE),
'response variable is missing from formula')
})

test_that("prior_only works", {
mod <- mvgam(y ~ s(season),
trend_model = AR(p = 2),
Expand Down Expand Up @@ -170,41 +277,6 @@ test_that("prior_only works", {
mod$model_file, fixed = TRUE)))
})

test_that("response variable must follow family-specific restrictions", {
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = gaus_data$data_train,
family = lognormal(),
run_model = FALSE),
'Values <= 0 not allowed for lognormal responses')

expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = gaus_data$data_train,
family = poisson(),
run_model = FALSE),
'Values < 0 not allowed for count family responses')
})

test_that("trend_model must be correctly specified", {
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR11',
data = beta_data$data_train,
family = betar(),
run_model = FALSE))
})

test_that("outcome variable must be present in data", {
data = data.frame(out = rnorm(100),
temp = rnorm(100),
time = 1:100)
expect_error(mod <- mvgam(formula = y ~ dynamic(temp, rho = 20),
data = data,
family = gaussian(),
run_model = FALSE),
'variable y not found in data')
})

test_that("time not required in data if this is a no trend model", {
data = data.frame(out = rnorm(100),
temp = rnorm(100))
Expand All @@ -215,39 +287,6 @@ test_that("time not required in data if this is a no trend model", {
expect_true(inherits(mod, 'mvgam_prefit'))
})

test_that("rho argument must be positive numeric", {
data = data.frame(out = rnorm(100),
temp = rnorm(100),
time = 1:100)
expect_error(mod <- mvgam(formula = out ~ dynamic(temp, rho = -1),
data = data,
family = gaussian(),
run_model = FALSE),
'Argument "rho" in dynamic() must be a positive value',
fixed = TRUE)
})

test_that("all series must have observations for all unique timepoints", {
data <- sim_mvgam()
data$data_train <- data$data_train[-2,]
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = data$data_train,
family = poisson(),
run_model = FALSE),
'One or more series in data is missing observations for one or more timepoints')

data <- sim_mvgam()
data$data_test <- data$data_test[-2,]
expect_error(mod <- mvgam(y ~ s(season),
trend_model = 'AR1',
data = data$data_train,
newdata = data$data_test,
family = poisson(),
run_model = FALSE),
'One or more series in newdata is missing observations for one or more timepoints')
})

test_that("median coefs should be stored in the mgcv object", {
expect_true(identical(unname(coef(mvgam:::mvgam_example2$mgcv_model)),
coef(mvgam:::mvgam_example2)[,2]))
Expand Down Expand Up @@ -277,36 +316,6 @@ test_that("empty obs formula allowed if trend_formula supplied", {
fixed = TRUE)))
})

test_that("missing values not allowed in predictors", {
# Include missing vals in training data
simdat <- sim_mvgam()
simdat$data_train$season[4] <- NA
expect_error(mvgam(y ~ s(season),
trend_model = 'GP',
data = simdat$data_train,
newdata = simdat$data_test,
run_model = FALSE))

# Include missing vals in testing data
simdat <- sim_mvgam()
simdat$data_test$season[4] <- NA
expect_error(mvgam(y ~ s(season),
data = simdat$data_train,
newdata = simdat$data_test,
run_model = FALSE))
})

test_that("series levels must match unique entries in series", {
# Include missing vals in training data
simdat <- sim_mvgam()
levels(simdat$data_train$series) <- paste0('series_', 1:6)
expect_error(mvgam(y ~ s(season),
trend_model = 'GP',
data = simdat$data_train,
newdata = simdat$data_test,
run_model = FALSE))
})

test_that("share_obs_params working", {
# Standard beta
mod <- mvgam(y ~ s(season, by = series),
Expand Down

0 comments on commit 9ced215

Please sign in to comment.