tests/testthat/test-compare_nested_models.R

comparison_female_male_cross_dyads <- dyads_cross |>
  dplyr::group_by(coupleID) |>
  dplyr::filter(dplyr::n_distinct(gender) == 2L) |>
  dplyr::ungroup()

test_that("comparison treats unregistered prefixed columns as original data", {
  data <- data.frame(
    dyad_id = c(1, 1, 2, 2),
    person_id = c("A", "B", "C", "D")
  )
  prepared <- prepare_dyad_data(
    data,
    dyad = dyad_id,
    member = person_id,
    seed = 123
  )
  changed <- prepared
  changed$.unregistered <- 1

  expect_error(
    validate_comparison_data(NULL, NULL, prepared, changed),
    "do not contain the same original columns"
  )
})

test_that("compare_nested_models compares reparameterized nested models", {
  skip_if_not_installed("glmmTMB")

  full_data <- prepare_dyad_data(
    dyads_cross,
    dyad = coupleID,
    member = personID,
    role = gender,
    seed = 123
  )
  restricted_data <- prepare_dyad_data(
    dyads_cross,
    dyad = coupleID,
    member = personID,
    role = gender,
    set_exchangeable_compositions = "male-female",
    pool_compositions = list(
      non_female_x_female = c("male-male", "female-male")
    ),
    seed = 123
  )

  full_model <- glmmTMB::glmmTMB(
    closeness ~ 0 +
      .is_female_x_male_female +
      .is_female_x_male_male +
      .is_female_x_female +
      .is_male_x_male +
      us(
        0 + .is_female_x_male_female + .is_female_x_male_male |
          coupleID
      ) +
      us(0 + .is_female_x_female | coupleID) +
      us(0 + .member_contrast_female_x_female_arbitrary | coupleID) +
      us(0 + .is_male_x_male | coupleID) +
      us(0 + .member_contrast_male_x_male_arbitrary | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = full_data
  )
  restricted_model <- glmmTMB::glmmTMB(
    closeness ~ 0 +
      .is_female_x_female +
      .is_non_female_x_female +
      us(0 + .is_female_x_female | coupleID) +
      us(0 + .member_contrast_female_x_female_arbitrary | coupleID) +
      us(0 + .is_non_female_x_female | coupleID) +
      us(0 + .member_contrast_non_female_x_female_arbitrary | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = restricted_data
  )

  comparison <- compare_nested_models(restricted_model, full_model)

  expect_s3_class(comparison, "anova")
  expect_equal(comparison$`Chi Df`[2], 5)
  expect_equal(
    comparison$Chisq[2],
    2 * as.numeric(logLik(full_model) - logLik(restricted_model))
  )
  expect_lt(comparison$`Pr(>Chisq)`[2], 0.001)
})

test_that("one prepared object supports both distinguishability models", {
  skip_if_not_installed("glmmTMB")

  comparison_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender,
    predictors = provided_support,
    include_arbitrary_member_contrast = TRUE,
    seed = 123
  )
  exchangeable_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender,
    predictors = provided_support,
    set_exchangeable_compositions = "female-male",
    seed = 123
  )

  full_model <- glmmTMB::glmmTMB(
    closeness ~ 0 +
      .is_female +
      .is_male +
      .is_female:.provided_support_actor +
      .is_male:.provided_support_actor +
      .is_female:.provided_support_partner +
      .is_male:.provided_support_partner +
      us(0 + .is_female + .is_male | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = comparison_data
  )
  restricted_model <- glmmTMB::glmmTMB(
    closeness ~
      .provided_support_actor +
      .provided_support_partner +
      us(1 | coupleID) +
      us(0 + .member_contrast_arbitrary | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = comparison_data
  )
  re_prepared_model <- update(
    restricted_model,
    data = exchangeable_data
  )

  comparison <- compare_nested_models(
    restricted_model,
    full_model
  )
  recovered <- recover_exchangeable_covariance(restricted_model)

  expect_s3_class(comparison, "anova")
  expect_equal(comparison$`Chi Df`[2], 4)
  expect_s3_class(recovered, "exchangeable_covariance")
  expect_equal(
    as.numeric(stats::logLik(restricted_model)),
    as.numeric(stats::logLik(re_prepared_model)),
    tolerance = 1e-8
  )
  expect_equal(
    glmmTMB::fixef(restricted_model),
    glmmTMB::fixef(re_prepared_model),
    tolerance = 1e-8
  )
  expect_equal(
    glmmTMB::VarCorr(restricted_model),
    glmmTMB::VarCorr(re_prepared_model),
    tolerance = 1e-8
  )
})

test_that("compare_nested_models requires exact original data", {
  skip_if_not_installed("glmmTMB")

  data_one <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender
  )
  changed <- comparison_female_male_cross_dyads
  changed$closeness[1] <- changed$closeness[1] + 1e-10
  data_two <- prepare_dyad_data(
    changed,
    dyad = coupleID,
    member = personID,
    role = gender
  )

  model_one <- glmmTMB::glmmTMB(
    closeness ~ 1 + us(1 | coupleID),
    data = data_one
  )
  model_two <- glmmTMB::glmmTMB(
    closeness ~ gender + us(1 | coupleID),
    data = data_two
  )

  expect_error(
    compare_nested_models(model_one, model_two),
    "Column `closeness` differs"
  )

  changed_metadata <- data_one
  attr(changed_metadata, "dyadMLM")$role <- "changed_role"
  expect_error(
    validate_comparison_data(
      model_one,
      model_one,
      data_one,
      changed_metadata
    ),
    "prepared datasets use different `role` variables"
  )
})

test_that("compare_nested_models checks source and fitted rows", {
  skip_if_not_installed("glmmTMB")

  complete_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender
  )
  first_dyad <- comparison_female_male_cross_dyads$coupleID[[1L]]
  fewer_rows <- comparison_female_male_cross_dyads[
    comparison_female_male_cross_dyads$coupleID != first_dyad,
  ]
  shorter_data <- prepare_dyad_data(
    fewer_rows,
    dyad = coupleID,
    member = personID,
    role = gender
  )
  predictor_missing <- comparison_female_male_cross_dyads
  predictor_missing$extra_predictor <- seq_len(nrow(predictor_missing))
  predictor_missing$extra_predictor[1] <- NA_real_
  row_data <- prepare_dyad_data(
    predictor_missing,
    dyad = coupleID,
    member = personID,
    role = gender
  )

  complete_restricted <- glmmTMB::glmmTMB(
    closeness ~ 1 + us(1 | coupleID),
    data = complete_data
  )
  shorter_full <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = shorter_data
  )
  row_restricted <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = row_data
  )
  row_full <- glmmTMB::glmmTMB(
    closeness ~ gender + extra_predictor,
    data = row_data
  )

  expect_error(
    compare_nested_models(complete_restricted, shorter_full),
    "different numbers of rows"
  )
  expect_error(
    compare_nested_models(row_restricted, row_full),
    "different observation rows"
  )
})

test_that("compare_nested_models supports and checks model families", {
  skip_if_not_installed("glmmTMB")

  data <- comparison_female_male_cross_dyads
  data$binary_outcome <- as.integer(
    data$closeness > stats::median(data$closeness, na.rm = TRUE)
  )
  restricted_data <- prepare_dyad_data(
    data,
    dyad = coupleID,
    member = personID,
    role = gender
  )
  full_data <- prepare_dyad_data(
    data,
    dyad = coupleID,
    member = personID,
    role = gender
  )

  restricted_model <- glmmTMB::glmmTMB(
    binary_outcome ~ 1,
    family = stats::binomial(),
    data = restricted_data
  )
  full_model <- glmmTMB::glmmTMB(
    binary_outcome ~ gender,
    family = stats::binomial(),
    data = full_data
  )
  gaussian_model <- glmmTMB::glmmTMB(
    binary_outcome ~ gender,
    family = stats::gaussian(),
    data = full_data
  )

  comparison <- compare_nested_models(restricted_model, full_model)

  expect_s3_class(comparison, "anova")
  expect_equal(comparison$`Chi Df`[2], 1)
  expect_equal(
    comparison$Chisq[2],
    2 * as.numeric(logLik(full_model) - logLik(restricted_model))
  )
  expect_error(
    compare_nested_models(restricted_model, gaussian_model),
    "same family and link"
  )
})

test_that("compare_nested_models compares APIM, DIM, and DSM models", {
  skip_if_not_installed("glmmTMB")

  distinguishable_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender,
    predictors = provided_support,
    model_types = c("apim", "dsm"),
    dsm_role_order = c("female", "male"),
    short_colnames = FALSE
  )
  exchangeable_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender,
    predictors = provided_support,
    model_types = c("apim", "dim"),
    set_exchangeable_compositions = "female-male",
    short_colnames = FALSE,
    seed = 123
  )

  distinguishable_dsm <- glmmTMB::glmmTMB(
    closeness ~ 1 +
      .provided_support_dyad_mean_gmc +
      .provided_support_within_dyad_diff +
      .dsm_role_contrast +
      .provided_support_dyad_mean_gmc:.dsm_role_contrast +
      .provided_support_within_dyad_diff:.dsm_role_contrast +
      us(1 + .dsm_role_contrast | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = distinguishable_data
  )
  exchangeable_apim <- glmmTMB::glmmTMB(
    closeness ~ 0 + .is_female_x_male +
      .provided_support_actor +
      .provided_support_partner +
      us(0 + .is_female_x_male | coupleID) +
      us(0 + .member_contrast_female_x_male_arbitrary | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = exchangeable_data
  )
  distinguishable_apim <- glmmTMB::glmmTMB(
    closeness ~ 0 +
      .is_female_x_male_female +
      .is_female_x_male_male +
      .is_female_x_male_female:.provided_support_actor +
      .is_female_x_male_male:.provided_support_actor +
      .is_female_x_male_female:.provided_support_partner +
      .is_female_x_male_male:.provided_support_partner +
      us(
        0 + .is_female_x_male_female + .is_female_x_male_male |
          coupleID
      ),
    dispformula = ~0,
    family = gaussian(),
    data = distinguishable_data
  )
  exchangeable_dim <- glmmTMB::glmmTMB(
    closeness ~ 1 +
      .provided_support_dyad_mean_gmc +
      .provided_support_within_dyad_dev +
      us(1 | coupleID) +
      us(0 + .member_contrast_female_x_male_arbitrary | coupleID),
    dispformula = ~0,
    family = gaussian(),
    data = exchangeable_data
  )

  dsm_comparison <- compare_nested_models(
    distinguishable_dsm,
    exchangeable_apim
  )
  apim_dim_comparison <- compare_nested_models(
    distinguishable_apim,
    exchangeable_dim
  )

  expect_equal(dsm_comparison$`Chi Df`[2], 4)
  expect_equal(apim_dim_comparison$`Chi Df`[2], 4)
  expect_equal(dsm_comparison$Chisq[2], apim_dim_comparison$Chisq[2])
})

test_that("compare_nested_models supports ordinary and mixed named data", {
  skip_if_not_installed("glmmTMB")

  plain_data_one <- comparison_female_male_cross_dyads
  plain_data_two <- plain_data_one
  prepared_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender
  )

  plain_smaller <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = plain_data_one
  )
  plain_larger <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = plain_data_two
  )
  prepared_larger <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = prepared_data
  )
  changed_data <- plain_data_two
  changed_data$provided_support[1] <- changed_data$provided_support[1] + 1
  changed_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = changed_data
  )
  inline_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = as.data.frame(plain_data_one)
  )

  plain_comparison <- compare_nested_models(plain_smaller, plain_larger)
  mixed_comparison <- compare_nested_models(plain_smaller, prepared_larger)

  expect_equal(plain_comparison$`Chi Df`[2], 1)
  expect_equal(
    plain_comparison$Chisq[2],
    2 * as.numeric(logLik(plain_larger) - logLik(plain_smaller))
  )
  expect_equal(
    mixed_comparison$Chisq[2],
    plain_comparison$Chisq[2]
  )
  expect_error(
    compare_nested_models(plain_smaller, changed_model),
    "Column `provided_support` differs"
  )
  expect_error(
    compare_nested_models(inline_model, plain_larger),
    "must have been fitted with a named data frame"
  )

  extra_column_data <- plain_data_two
  extra_column_data$unused_column <- 1
  expect_error(
    validate_comparison_data(
      plain_smaller,
      plain_larger,
      plain_data_one,
      extra_column_data
    ),
    "do not contain the same original columns"
  )
})

test_that("compare_nested_models sorts models and agrees with anova.glmmTMB", {
  skip_if_not_installed("glmmTMB")

  model_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender
  )
  restricted_model <- glmmTMB::glmmTMB(
    closeness ~ 1 + us(1 | coupleID),
    data = model_data
  )
  full_model <- glmmTMB::glmmTMB(
    closeness ~ gender + us(1 | coupleID),
    data = model_data
  )

  reference <- stats::anova(restricted_model, full_model)
  comparison <- compare_nested_models(restricted_model, full_model)
  reversed <- compare_nested_models(full_model, restricted_model)

  expect_equal(comparison$Chisq[2], reference$Chisq[2])
  expect_equal(comparison$`Chi Df`[2], reference$`Chi Df`[2])
  expect_equal(comparison$`Pr(>Chisq)`[2], reference$`Pr(>Chisq)`[2])
  expect_identical(reversed, comparison)
  expect_identical(row.names(comparison), c("restricted_model", "full_model"))
  expect_identical(
    comparison$`Pr(>Chisq)`[2] < 0.05,
    reference$`Pr(>Chisq)`[2] < 0.05
  )
  expect_error(
    compare_nested_models(restricted_model, restricted_model),
    "different numbers of estimated parameters"
  )

  printed <- capture.output(print(comparison))
  note <- match(
    paste0(
      "Assumes mathematical nesting and an appropriate chi-squared ",
      "reference distribution."
    ),
    printed
  )
  expect_equal(printed[note + 1L], "")
  expect_match(
    tail(printed, 1L),
    "`full_model` fits better than `restricted_model`",
    fixed = TRUE
  )
  formatted_p <- format.pval(
    comparison$`Pr(>Chisq)`[2],
    digits = 3,
    eps = 0.001
  )
  if (startsWith(formatted_p, "<")) {
    formatted_p <- paste("p <", substring(formatted_p, 2L))
  } else {
    formatted_p <- paste("p =", formatted_p)
  }
  expect_match(
    tail(printed, 1L),
    formatted_p,
    fixed = TRUE
  )

  set.seed(195)
  model_data$noise <- stats::rnorm(nrow(model_data))
  smaller_model <- glmmTMB::glmmTMB(closeness ~ 1, data = model_data)
  larger_model <- glmmTMB::glmmTMB(closeness ~ noise, data = model_data)
  no_clear_improvement <- compare_nested_models(smaller_model, larger_model)
  printed <- capture.output(print(no_clear_improvement))

  expect_gt(no_clear_improvement$`Pr(>Chisq)`[2], 0.05)
  expect_match(
    tail(printed, 1L),
    "finds no clear improvement from `smaller_model` to `larger_model`",
    fixed = TRUE
  )
  expect_match(
    tail(printed, 1L),
    "This does not establish equal fit",
    fixed = TRUE
  )
  expect_false(
    grepl(
      "prefer `smaller_model` for parsimony",
      tail(printed, 1L),
      fixed = TRUE
    )
  )
  expect_match(
    tail(printed, 1L),
    "Conclusion (5% level)",
    fixed = TRUE
  )
  expect_match(
    tail(printed, 1L),
    format.pval(
      no_clear_improvement$`Pr(>Chisq)`[2],
      digits = 3,
      eps = 0.001
    ),
    fixed = TRUE
  )

  p_value <- no_clear_improvement$`Pr(>Chisq)`[2]
  higher_alpha <- (p_value + 1) / 2
  comparison_at_higher_alpha <- compare_nested_models(
    smaller_model,
    larger_model,
    alpha = higher_alpha
  )

  expect_equal(attr(comparison_at_higher_alpha, "alpha"), higher_alpha)
  expect_match(
    tail(capture.output(print(comparison_at_higher_alpha)), 1L),
    "provides evidence that `larger_model` fits better than `smaller_model`",
    fixed = TRUE
  )
})

test_that("compare_nested_models validates alpha", {
  invalid_alpha_values <- list(
    NULL,
    NA_real_,
    Inf,
    0,
    1,
    c(0.01, 0.05),
    "0.05"
  )

  for (invalid_alpha in invalid_alpha_values) {
    expect_error(
      compare_nested_models(NULL, NULL, alpha = invalid_alpha),
      "`alpha` must be one finite number between 0 and 1",
      fixed = TRUE
    )
  }
})

test_that("compare_nested_models supports fits without standard errors", {
  skip_if_not_installed("glmmTMB")

  model_data <- comparison_female_male_cross_dyads
  restricted_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = model_data,
    se = FALSE
  )
  full_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = model_data,
    se = FALSE
  )

  comparison <- compare_nested_models(restricted_model, full_model)

  expect_null(restricted_model$sdr$pdHess)
  expect_null(full_model$sdr$pdHess)
  expect_s3_class(comparison, "anova")
  expect_equal(comparison$`Chi Df`[2L], 1)
  expect_true(is.finite(comparison$Chisq[2L]))
})

test_that("compare_nested_models rejects transformed and changed outcomes", {
  skip_if_not_installed("glmmTMB")

  model_data <- prepare_dyad_data(
    comparison_female_male_cross_dyads,
    dyad = coupleID,
    member = personID,
    role = gender
  )
  restricted_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = model_data
  )
  transformed_model <- glmmTMB::glmmTMB(
    I(closeness) ~ gender,
    data = model_data
  )

  model_data$alternate_outcome <- model_data$closeness
  alternate_outcome_model <- glmmTMB::glmmTMB(
    alternate_outcome ~ gender,
    data = model_data
  )

  expect_error(
    compare_nested_models(restricted_model, transformed_model),
    "Only an untransformed outcome"
  )
  expect_error(
    compare_nested_models(restricted_model, alternate_outcome_model),
    "use different outcome variables"
  )

  model_data$closeness <- model_data$closeness + 1
  changed_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = model_data
  )
  expect_error(
    compare_nested_models(restricted_model, changed_model),
    "fitted outcome values differ"
  )
})

test_that("compare_nested_models recovers local model data", {
  skip_if_not_installed("glmmTMB")

  fit_models <- function() {
    local_data <- prepare_dyad_data(
      comparison_female_male_cross_dyads,
      dyad = coupleID,
      member = personID,
      role = gender
    )
    list(
      restricted = glmmTMB::glmmTMB(closeness ~ 1, data = local_data),
      full = glmmTMB::glmmTMB(closeness ~ gender, data = local_data)
    )
  }

  models <- fit_models()
  comparison <- compare_nested_models(models$restricted, models$full)

  expect_s3_class(comparison, "anova")
})

test_that("compare_nested_models rejects unsuitable fitted models", {
  skip_if_not_installed("glmmTMB")

  model_data <- comparison_female_male_cross_dyads
  valid_model <- glmmTMB::glmmTMB(closeness ~ 1, data = model_data)

  expect_error(
    compare_nested_models(stats::lm(closeness ~ 1, data = model_data), valid_model),
    "`model1` must be a fitted `glmmTMB` model"
  )

  reml_model <- valid_model
  reml_model$modelInfo$REML <- TRUE
  expect_error(
    compare_nested_models(valid_model, reml_model),
    "`model2` model was fitted with REML"
  )

  unconverged_model <- valid_model
  unconverged_model$fit$convergence <- 1L
  expect_error(
    compare_nested_models(valid_model, unconverged_model),
    "`model2` model did not converge"
  )

  indefinite_hessian_model <- valid_model
  indefinite_hessian_model$sdr$pdHess <- FALSE
  expect_error(
    compare_nested_models(valid_model, indefinite_hessian_model),
    "`model2` model has a non-positive-definite Hessian matrix"
  )

  unavailable_data_model <- valid_model
  unavailable_data_model$call$data <- as.name("unavailable_model_data")
  expect_error(
    comparison_model_data(
      unavailable_data_model,
      parent.frame(),
      "model1"
    ),
    "Could not recover the named data frame used by `model1`"
  )
})

test_that("compare_nested_models rejects models fitted with priors", {
  skip_if_not_installed("glmmTMB")

  model_data <- comparison_female_male_cross_dyads
  intercept_prior <- function(prior) {
    data.frame(prior = prior, class = "fixef", coef = "(Intercept)")
  }
  restricted_model <- glmmTMB::glmmTMB(closeness ~ 1, data = model_data)
  full_model <- glmmTMB::glmmTMB(closeness ~ gender, data = model_data)
  prior_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = model_data,
    priors = intercept_prior("normal(0,1)")
  )
  restricted_prior_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = model_data,
    priors = intercept_prior("normal(0,1)")
  )
  restricted_other_prior_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    data = model_data,
    priors = intercept_prior("normal(0,10)")
  )

  expect_error(
    compare_nested_models(restricted_model, prior_model),
    "`model2` model was fitted with priors"
  )
  expect_error(
    compare_nested_models(prior_model, restricted_model),
    "`model1` model was fitted with priors"
  )
  # Identical priors still give a maximum a posteriori fit.
  expect_error(
    compare_nested_models(restricted_prior_model, prior_model),
    "`model1` model was fitted with priors"
  )
  expect_error(
    compare_nested_models(restricted_other_prior_model, prior_model),
    "`model1` model was fitted with priors"
  )

  # NULL and empty prior specifications fit without priors.
  null_prior_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = model_data,
    priors = NULL
  )
  empty_prior_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = model_data,
    priors = data.frame(
      prior = character(),
      class = character(),
      coef = character()
    )
  )
  expected <- compare_nested_models(restricted_model, full_model)

  expect_equal(
    compare_nested_models(restricted_model, null_prior_model)$Chisq,
    expected$Chisq
  )
  expect_equal(
    compare_nested_models(restricted_model, empty_prior_model)$Chisq,
    expected$Chisq
  )
})

test_that("compare_nested_models checks weights and offsets", {
  skip_if_not_installed("glmmTMB")

  data <- comparison_female_male_cross_dyads
  data$observation_weight <- 2
  data$observation_offset <- 0.25
  model_data <- prepare_dyad_data(
    data,
    dyad = coupleID,
    member = personID,
    role = gender
  )

  full_model <- glmmTMB::glmmTMB(
    closeness ~ gender,
    data = model_data
  )
  weighted_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    weights = observation_weight,
    data = model_data
  )
  offset_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    offset = observation_offset,
    data = model_data
  )
  zero_inflation_offset_model <- glmmTMB::glmmTMB(
    closeness ~ 1,
    ziformula = ~0 + offset(observation_offset),
    data = model_data
  )

  expect_error(
    compare_nested_models(full_model, weighted_model),
    "different observation weights"
  )
  expect_error(
    compare_nested_models(full_model, offset_model),
    "different offsets"
  )
  expect_error(
    compare_nested_models(offset_model, zero_inflation_offset_model),
    "different offsets"
  )
})

Try the dyadMLM package in your browser

Any scripts or data that you put into this service are public.

dyadMLM documentation built on Oct. 28, 2026, 1:06 a.m.