R/train_GLM.R

In rtemis: Schema-Driven Machine Learning and Visualization

# train_GLM.R
# ::rtemis::
# 2025- EDG rtemis.org

# %% train_.GLMHyperparameters ----
#' Train a GLM model
#'
#' Train a GLM model using `stats::glm`.
#'
#' @details
#' `stats::glm` does not work in the presence of missing values.
#' This function uses the formula interface to `glm` to train a GLM model.
#' No preprocessing is needed.
#'
#' @param x tabular data: Training set.
#' @param weights Numeric vector: Case weights.
#' @param hyperparameters `GLMHyperparameters` object: make using [setup_GLM].
#' @param verbosity Integer: If > 0, print messages.
#'
#' @return GLM model.
#'
#' @author EDG
#' @keywords internal
#' @noRd
method(train_, GLMHyperparameters) <- function(
  hyperparameters,
  x,
  weights = NULL,
  dat_validation = NULL,
  execution_config = setup_FutureExecution(),
  verbosity = 1L
) {
  # Data ----
  check_supervised(
    x = x,
    allow_missing = FALSE,
    verbosity = verbosity
  )

  if (is.null(weights)) {
    weights <- rep(1, NROW(x))
  }

  type <- supervised_type(x)
  if (type == "Classification") {
    n_classes <- nlevels(outcome(x))
    if (n_classes > 2L) {
      rtemis.core::abort(
        "GLM does not support multiclass classification.",
        class = "rtemis_unsupported_error"
      )
    }
  } else {
    n_classes <- NA_integer_
  }

  # Formula ----
  formula <- as.formula(
    paste(
      names(x)[ncol(x)],
      "~ ."
    )
  )

  # Train ----
  family <- if (type == "Regression") {
    gaussian()
  } else if (type == "Classification") {
    binomial()
  }
  model <- glm(
    formula = formula,
    family = family,
    data = x,
    weights = weights
  )
  check_inherits(model, "glm")
  list(model = model, preprocessor = NULL)
} # /rtemis::train_.GLMHyperparameters


# %% predict_super.class_glm ----
#' Predict from GLM model
#'
#' @param model GLM model.
#' @param newdata data.frame or similar: Data to predict on.
#'
#' @keywords internal
#' @noRd
method(predict_super, class_glm) <- function(
  model,
  newdata,
  type = NULL,
  execution_config = NULL,
  verbosity = 0L
) {
  predict(model, newdata = newdata, type = "response")
} # /rtemis::predict_super.glm


# %% varimp_super.class_glm ----
#' Get coefficients from GLM model
#'
#' @param model GLM model.
#'
#' @keywords internal
#' @noRd
method(varimp_super, class_glm) <- function(model) {
  .coef <- coef(model)
  VariableImportance(
    measures = list(
      Coefficient = importance_measure(
        names(.coef),
        unname(.coef),
        kind = "coefficient",
        signed = TRUE,
        scale_dependent = TRUE,
        direction = "absolute",
        description = paste0(
          "Estimated coefficient of each model term, including the intercept, ",
          "on the scale of the linear predictor, per unit of the term as the ",
          "model received it."
        )
      )
    )
  )
} # /rtemis::varimp_super.glm


# %% se_super.class_glm ----
#' Get Standard Errors from GLM model
#'
#' @param model GLM model.
#' @param newdata data.frame or similar: Data to predict on.
#'
#' @author EDG
#' @keywords internal
#' @noRd
method(se_super, class_glm) <- function(model, newdata) {
  predict(model, newdata = newdata, se.fit = TRUE)[["se.fit"]]
} # /rtemis::se_super.glm


# %% explain_super.class_glm ----
#' LinearSHAP contributions from a GLM
#'
#' Exact and closed-form: the model is linear on its link scale, so
#' `phi_j = beta_j * (x_j - E[x_j])`. For a binomial family that scale is the
#' logit, which is why the contributions do not sum to the fitted probability.
#'
#' @param model `glm` object.
#' @param newdata tabular data: Cases to explain.
#' @param background tabular data: Reference cases.
#' @param estimator Character: Resolved estimator.
#' @param perturbation Character: Resolved value function.
#' @param scale Character: Scale the contributions are additive on.
#' @param type Character: "Regression" or "Classification".
#' @param verbosity Integer: Verbosity level.
#'
#' @return List with `phi`, `baseline`, `predicted` and `exact`.
#'
#' @keywords internal
#' @noRd
method(explain_super, class_glm) <- function(
  model,
  newdata,
  background,
  estimator,
  perturbation,
  scale,
  type,
  verbosity = 0L
) {
  check_shap_linear(estimator, perturbation, "GLM")
  shap_require_background(background, "LinearSHAP")
  model_matrix_shap(
    newdata = newdata,
    background = background,
    formula_terms = stats::delete.response(stats::terms(model)),
    coefficients = matrix(stats::coef(model), ncol = 1L),
    margin_fn = function(design) {
      matrix(
        predict(model, newdata = as.data.frame(newdata), type = "link"),
        ncol = 1L
      )
    },
    xlev = model[["xlevels"]],
    contrasts.arg = model[["contrasts"]],
    label = "LinearSHAP"
  )
} # /rtemis::explain_super.glm