From 6efcf3ee3ca46b45cf4be6645543b5defe8c6b36 Mon Sep 17 00:00:00 2001
From: Emil Hvitfeldt <emilhhvitfeldt@gmail.com>
Date: Sun, 23 Aug 2026 21:43:50 -0700
Subject: [PATCH] Add generics describing what fitted expressions compute

tidypredict_fit() returns expressions but says nothing about what they
mean, and the meaning is not recoverable from their shape. Two cases
where getting it wrong is silent:

  A LiblineaR logistic regression and a LiblineaR SVM classifier both
  return one expression. The first is a probability in [0, 1]; the second
  is an uncalibrated decision value that straddles 0. Cutting the second
  at 0.5 as though it were a probability misclassifies every row whose
  value lies between 0 and 0.5.

  A multiclass probability list and a quantreg::rq() fit with several tau
  are both named lists of expressions of the same length. In the first
  the values sum to one; in the second they are unrelated predictions.

Add tidypredict_output_type(), tidypredict_outcome_levels() and
tidypredict_normalized(), asked of the model rather than of the result.
Attributes on the result were the obvious alternative and do not work:
they are dropped by the subsetting, lapply() and unlist() that callers
apply to a multiclass list, and tidypredict_save() serialises to YAML so
they do not survive a round trip either.

Only four output types exist. There is no "logit": every multiclass list
goes through one shared expr_softmax(), so tidypredict_normalized() is
always TRUE or NA and never FALSE. It is worth keeping so callers can
rely on that rather than having to know it.

Two cases are less obvious than they look. rpart, ctree, C5.0 and bagger
classification return a class label and no probability at all, so they
are "class" rather than "prob". xgboost's binary:hinge objective is also
"class": it is wrapped in as.numeric(score >= 0) and so takes only the
values 0 and 1, which makes using it as a numeric prediction wrong even
though its type is numeric.

There are 44 tidypredict_fit() methods but only 10 parsed-model types, so
the .default methods route through parse_model() the way
tidypredict_fit.default() does, and most methods sit on the pm_* classes.
A backend whose parsed form does not record what is needed gets a method
on the fitted class: rpart and party do not record their mode, ksvm and
glm do not record their response levels once parsed, and parse_model.rqs()
returns a bare list rather than a parsed model.
---
 NAMESPACE                                |  66 +++
 NEWS.md                                  |   6 +
 R/model-C5.0.R                           |  22 +
 R/model-bagger.R                         |  30 ++
 R/model-bart.R                           |  21 +
 R/model-catboost.R                       |  40 ++
 R/model-glm.R                            |  20 +
 R/model-glmnet.R                         |  11 +
 R/model-h2o.R                            |  63 +++
 R/model-ksvm.R                           |  24 +-
 R/model-liblinear.R                      |  56 +++
 R/model-lightgbm.R                       |  39 ++
 R/model-nnet.R                           |  32 ++
 R/model-partykit.R                       |  36 ++
 R/model-rpart.R                          |  33 ++
 R/model-rq.R                             |  26 +
 R/model-xgboost.R                        |  41 ++
 R/output-metadata.R                      | 343 ++++++++++++++
 man/tidypredict_metadata.Rd              |  99 ++++
 tests/testthat/_snaps/output-metadata.md |  16 +
 tests/testthat/test-output-metadata.R    | 577 +++++++++++++++++++++++
 21 files changed, 1600 insertions(+), 1 deletion(-)
 create mode 100644 R/output-metadata.R
 create mode 100644 man/tidypredict_metadata.Rd
 create mode 100644 tests/testthat/_snaps/output-metadata.md
 create mode 100644 tests/testthat/test-output-metadata.R

diff --git a/NAMESPACE b/NAMESPACE
index 71c863ce..b441dfcc 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -125,6 +125,69 @@ S3method(tidypredict_n_trees,lgb.Booster)
 S3method(tidypredict_n_trees,randomForest)
 S3method(tidypredict_n_trees,ranger)
 S3method(tidypredict_n_trees,xgb.Booster)
+S3method(tidypredict_normalized,C5.0)
+S3method(tidypredict_normalized,H2OBinomialModel)
+S3method(tidypredict_normalized,H2OMultinomialModel)
+S3method(tidypredict_normalized,H2ORegressionModel)
+S3method(tidypredict_normalized,LiblineaR)
+S3method(tidypredict_normalized,bagger)
+S3method(tidypredict_normalized,default)
+S3method(tidypredict_normalized,party)
+S3method(tidypredict_normalized,pm_bart)
+S3method(tidypredict_normalized,pm_catboost)
+S3method(tidypredict_normalized,pm_lgb)
+S3method(tidypredict_normalized,pm_multiclass_regression)
+S3method(tidypredict_normalized,pm_naive_bayes)
+S3method(tidypredict_normalized,pm_nnet)
+S3method(tidypredict_normalized,pm_nullmodel_classification)
+S3method(tidypredict_normalized,pm_regression)
+S3method(tidypredict_normalized,pm_tree)
+S3method(tidypredict_normalized,pm_xgb)
+S3method(tidypredict_normalized,rpart)
+S3method(tidypredict_normalized,rqs)
+S3method(tidypredict_outcome_levels,C5.0)
+S3method(tidypredict_outcome_levels,H2OBinomialModel)
+S3method(tidypredict_outcome_levels,H2OMultinomialModel)
+S3method(tidypredict_outcome_levels,H2ORegressionModel)
+S3method(tidypredict_outcome_levels,LiblineaR)
+S3method(tidypredict_outcome_levels,bagger)
+S3method(tidypredict_outcome_levels,default)
+S3method(tidypredict_outcome_levels,glm)
+S3method(tidypredict_outcome_levels,ksvm)
+S3method(tidypredict_outcome_levels,lognet)
+S3method(tidypredict_outcome_levels,party)
+S3method(tidypredict_outcome_levels,pm_bart)
+S3method(tidypredict_outcome_levels,pm_catboost)
+S3method(tidypredict_outcome_levels,pm_lgb)
+S3method(tidypredict_outcome_levels,pm_multiclass_regression)
+S3method(tidypredict_outcome_levels,pm_naive_bayes)
+S3method(tidypredict_outcome_levels,pm_nnet)
+S3method(tidypredict_outcome_levels,pm_nullmodel_classification)
+S3method(tidypredict_outcome_levels,pm_regression)
+S3method(tidypredict_outcome_levels,pm_tree)
+S3method(tidypredict_outcome_levels,pm_xgb)
+S3method(tidypredict_outcome_levels,rpart)
+S3method(tidypredict_outcome_levels,rqs)
+S3method(tidypredict_output_type,C5.0)
+S3method(tidypredict_output_type,H2OBinomialModel)
+S3method(tidypredict_output_type,H2OMultinomialModel)
+S3method(tidypredict_output_type,H2ORegressionModel)
+S3method(tidypredict_output_type,LiblineaR)
+S3method(tidypredict_output_type,bagger)
+S3method(tidypredict_output_type,default)
+S3method(tidypredict_output_type,party)
+S3method(tidypredict_output_type,pm_bart)
+S3method(tidypredict_output_type,pm_catboost)
+S3method(tidypredict_output_type,pm_lgb)
+S3method(tidypredict_output_type,pm_multiclass_regression)
+S3method(tidypredict_output_type,pm_naive_bayes)
+S3method(tidypredict_output_type,pm_nnet)
+S3method(tidypredict_output_type,pm_nullmodel_classification)
+S3method(tidypredict_output_type,pm_regression)
+S3method(tidypredict_output_type,pm_tree)
+S3method(tidypredict_output_type,pm_xgb)
+S3method(tidypredict_output_type,rpart)
+S3method(tidypredict_output_type,rqs)
 S3method(tidypredict_test,"_xgb.Booster")
 S3method(tidypredict_test,C5.0)
 S3method(tidypredict_test,H2OBinomialModel)
@@ -201,6 +264,9 @@ export(tidypredict_fit)
 export(tidypredict_interval)
 export(tidypredict_load)
 export(tidypredict_n_trees)
+export(tidypredict_normalized)
+export(tidypredict_outcome_levels)
+export(tidypredict_output_type)
 export(tidypredict_save)
 export(tidypredict_sql)
 export(tidypredict_sql_interval)
diff --git a/NEWS.md b/NEWS.md
index d23f1592..dc94b96f 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,7 +1,13 @@
 # tidypredict (development version)
 
+- New generics describe what a model's fitted expressions compute, which the expressions themselves do not say: `tidypredict_output_type()` returns one of `"numeric"`, `"prob"`, `"decision"` or `"class"`, `tidypredict_outcome_levels()` returns the outcome levels in model order, and `tidypredict_normalized()` reports whether per-level probabilities already sum to one. See `?tidypredict_metadata`. (#435)
+
+- The distinctions these generics record cannot be recovered from the shape of a `tidypredict_fit()` result. A `LiblineaR` SVM classifier and a `LiblineaR` logistic regression both return a single expression, but the first is an uncalibrated decision value whose sign picks the class, so thresholding it at 0.5 as though it were a probability is wrong. A multiclass probability list and a `quantreg::rq()` fit with several `tau` are both named lists of expressions of the same length. (#435)
+
 - Fixed `cforest()` models failing with "'language' object cannot be coerced to type 'symbol'" under partykit 1.3-0. That release added a shim to partykit's methods that identifies the caller with `as.name()`, which errors when the generic is reached as `partykit::gettree()`. (#434)
 
+- `tidypredict_output_type()`, `tidypredict_outcome_levels()` and `tidypredict_normalized()` describe what a model's fitted expressions actually compute: a number, a probability, an uncalibrated decision value or a class label, which outcome levels they are named for, and whether a per-level list already sums to one. None of it is recoverable from the shape of the result, so a caller generating code from `tidypredict_fit()` no longer has to keep its own list of which backend produces which shape. See `?tidypredict_metadata`. (#433)
+
 - New generics expose the pieces `tidypredict_fit()` is assembled from, so that packages generating their own code from a fitted model can reuse tidypredict's parsing: `tidypredict_trees()` returns per-tree expressions, `tidypredict_class_trees()` returns per-tree expressions for each outcome level, `tidypredict_class_exprs()` returns one finished expression per outcome level, and `tidypredict_n_trees()` returns the number of trees. See `?tidypredict_extractors`. (#433)
 
 - The eleven `.extract_*()` functions are deprecated in favour of those generics. They were exported but documented as internal, and each is now a thin wrapper that warns. Two of them change return type under the new names: `.extract_earth_multiclass()` and `.extract_glmnet_multiclass()` returned deparsed strings, while `tidypredict_class_exprs()` returns language objects like every other extractor. (#433)
diff --git a/R/model-C5.0.R b/R/model-C5.0.R
index ebeabbee..d0069754 100644
--- a/R/model-C5.0.R
+++ b/R/model-C5.0.R
@@ -1092,3 +1092,25 @@ build_tree_formula.pm_tree_C5.0 <- function(model) {
   }
   c50_with_na_descent(build_tree_formula_single(model), model$tree_info)
 }
+
+# Output metadata ---------------------------------
+
+# C5.0 is classification only, and the fit votes the trees into a single class
+# label rather than a probability.
+#' @export
+tidypredict_output_type.C5.0 <- function(x, ...) {
+  rlang::check_dots_empty()
+  "class"
+}
+
+#' @export
+tidypredict_outcome_levels.C5.0 <- function(x, ...) {
+  rlang::check_dots_empty()
+  x$levels
+}
+
+#' @export
+tidypredict_normalized.C5.0 <- function(x, ...) {
+  rlang::check_dots_empty()
+  NA
+}
diff --git a/R/model-bagger.R b/R/model-bagger.R
index 4a050ad8..90edb663 100644
--- a/R/model-bagger.R
+++ b/R/model-bagger.R
@@ -197,3 +197,33 @@ bagger_check_model <- function(model, call = rlang::caller_env()) {
 build_tree_formula.pm_tree_bagger <- function(model) {
   bagger_build_formula(model)
 }
+
+# Output metadata ---------------------------------
+
+# `bagger_classes()` is the same mode signal `parse_model.bagger()` uses:
+# `NULL` for regression, the outcome levels for classification. A classification
+# ensemble averages class probabilities and then picks the largest, so the fit
+# is a class label rather than a probability.
+#' @export
+tidypredict_output_type.bagger <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (is.null(bagger_classes(x))) {
+    return("numeric")
+  }
+  "class"
+}
+
+#' @export
+tidypredict_outcome_levels.bagger <- function(x, ...) {
+  rlang::check_dots_empty()
+  bagger_classes(x)
+}
+
+#' @export
+tidypredict_normalized.bagger <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # A single expression either way, so there are no per-level values to sum.
+  NA
+}
diff --git a/R/model-bart.R b/R/model-bart.R
index 9b884cb3..b80254f2 100644
--- a/R/model-bart.R
+++ b/R/model-bart.R
@@ -257,3 +257,24 @@ tidypredict_test.bart <- function(
     y_scale = diff(range(model$y))
   )
 }
+
+# Output metadata ---------------------------------
+
+# {dbarts} models are supported for regression only.
+#' @export
+tidypredict_output_type.pm_bart <- function(x, ...) {
+  rlang::check_dots_empty()
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.pm_bart <- function(x, ...) {
+  rlang::check_dots_empty()
+  NULL
+}
+
+#' @export
+tidypredict_normalized.pm_bart <- function(x, ...) {
+  rlang::check_dots_empty()
+  NA
+}
diff --git a/R/model-catboost.R b/R/model-catboost.R
index 321a3c3f..e99c3fc8 100644
--- a/R/model-catboost.R
+++ b/R/model-catboost.R
@@ -919,3 +919,43 @@ tidypredict_n_trees.catboost.Model <- function(x, ...) {
   # round, so `niter` is the number of rounds rather than the number of trees.
   length(tidypredict_trees(x))
 }
+
+# Output metadata ---------------------------------
+
+# The same objective groups `build_fit_formula_catboost_nested()` switches on.
+catboost_parsed_objective <- function(x) {
+  x$general$params$objective %||% "RMSE"
+}
+
+#' @export
+tidypredict_output_type.pm_catboost <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  objective <- catboost_parsed_objective(x)
+  if (
+    objective %in%
+      c(catboost_multiclass_objectives, catboost_sigmoid_objectives)
+  ) {
+    return("prob")
+  }
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.pm_catboost <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # CatBoost is fit on integer labels, and the multiclass expressions come back
+  # named `class_0`, `class_1` and so on, which are positions, not levels.
+  NULL
+}
+
+#' @export
+tidypredict_normalized.pm_catboost <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (catboost_parsed_objective(x) %in% catboost_multiclass_objectives) {
+    return(TRUE)
+  }
+  NA
+}
diff --git a/R/model-glm.R b/R/model-glm.R
index 2807a824..fe631ad7 100644
--- a/R/model-glm.R
+++ b/R/model-glm.R
@@ -29,3 +29,23 @@ te_interval_glm <- function(parsedmodel, interval = 0.95) {
   }
   intervals
 }
+
+# Output metadata ---------------------------------
+
+# The parsed form is coefficients only, so the response levels have to come off
+# the fitted object. A binomial glm fit on a 0/1 numeric records none, which is
+# the "did not retain the levels" case.
+#' @export
+tidypredict_outcome_levels.glm <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (!identical(x$family$family, "binomial")) {
+    return(NULL)
+  }
+
+  response <- stats::model.frame(x)[[1]]
+  if (!is.factor(response)) {
+    return(NULL)
+  }
+  levels(response)
+}
diff --git a/R/model-glmnet.R b/R/model-glmnet.R
index 9a18da4c..d60c0e0f 100644
--- a/R/model-glmnet.R
+++ b/R/model-glmnet.R
@@ -210,3 +210,14 @@ tidypredict_class_exprs.multnet <- function(x, ..., penalty = NULL) {
   names(eqs) <- class_names
   eqs
 }
+
+# Output metadata ---------------------------------
+
+# A `multnet` fit parses to a multiclass parsed model, which carries its own
+# levels. A binary `lognet` parses to a single logistic expression, so its
+# levels only exist on the fitted object.
+#' @export
+tidypredict_outcome_levels.lognet <- function(x, ...) {
+  rlang::check_dots_empty()
+  as.character(x$classnames)
+}
diff --git a/R/model-h2o.R b/R/model-h2o.R
index 4c66d3ca..c8187b01 100644
--- a/R/model-h2o.R
+++ b/R/model-h2o.R
@@ -208,3 +208,66 @@ tidypredict_test.H2OMultinomialModel <- function(
 h2o_test_results <- function(base, te, threshold) {
   test_results_numeric(base, te, threshold)
 }
+
+# Output metadata ---------------------------------
+
+# H2O models have no `parse_model()` method, so the metadata has to be read off
+# the model handle. The three classes are exactly the three modes.
+#' @export
+tidypredict_output_type.H2ORegressionModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.H2ORegressionModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  NULL
+}
+
+#' @export
+tidypredict_normalized.H2ORegressionModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  NA
+}
+
+#' @export
+tidypredict_output_type.H2OBinomialModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  "prob"
+}
+
+#' @export
+tidypredict_outcome_levels.H2OBinomialModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  h2o_response_domain(x)
+}
+
+#' @export
+tidypredict_normalized.H2OBinomialModel <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # `tidypredict_fit()` returns the probability of the second domain level
+  # only, so there are no per-level values to sum.
+  NA
+}
+
+#' @export
+tidypredict_output_type.H2OMultinomialModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  "prob"
+}
+
+#' @export
+tidypredict_outcome_levels.H2OMultinomialModel <- function(x, ...) {
+  rlang::check_dots_empty()
+  h2o_response_domain(x)
+}
+
+#' @export
+tidypredict_normalized.H2OMultinomialModel <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # `expr_softmax()` over the domain.
+  TRUE
+}
diff --git a/R/model-ksvm.R b/R/model-ksvm.R
index bf14f6b7..d0faea6a 100644
--- a/R/model-ksvm.R
+++ b/R/model-ksvm.R
@@ -5,6 +5,10 @@
 #' @export
 parse_model.ksvm <- function(model) parse_model_ksvm(model)
 
+ksvm_is_classification <- function(model) {
+  kernlab::type(model) %in% c("C-svc", "nu-svc", "C-bsvc")
+}
+
 parse_model_ksvm <- function(model, call = rlang::caller_env()) {
   acceptable_formula(model)
 
@@ -19,7 +23,7 @@ parse_model_ksvm <- function(model, call = rlang::caller_env()) {
   }
 
   svm_type <- kernlab::type(model)
-  is_classification <- svm_type %in% c("C-svc", "nu-svc", "C-bsvc")
+  is_classification <- ksvm_is_classification(model)
   is_regression <- svm_type %in% c("eps-svr", "nu-svr", "eps-bsvr")
   if (!is_classification && !is_regression) {
     cli::cli_abort(
@@ -287,3 +291,21 @@ tidypredict_test.ksvm <- function(
     model@kcall
   )
 }
+
+# Output metadata ---------------------------------
+
+# `parse_model_ksvm()` folds Platt scaling into the glm logit machinery, so
+# `tidypredict_output_type.pm_regression()` already reports "prob" for a
+# classifier and "numeric" for an SVR. Only the levels need the fitted object.
+#' @export
+tidypredict_outcome_levels.ksvm <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # `lev()` holds the sorted response values for a regression, so read the mode
+  # off `type()` the way `parse_model_ksvm()` does. A classifier is binary
+  # because `parse_model_ksvm()` rejects anything else.
+  if (!ksvm_is_classification(x)) {
+    return(NULL)
+  }
+  as.character(kernlab::lev(x))
+}
diff --git a/R/model-liblinear.R b/R/model-liblinear.R
index 3ed93879..53f898ff 100644
--- a/R/model-liblinear.R
+++ b/R/model-liblinear.R
@@ -95,6 +95,62 @@ parse_model_liblinear <- function(model, call = rlang::caller_env()) {
   as_parsed_model(pm)
 }
 
+# Output metadata ---------------------------------
+
+# This is the model class the metadata generics exist for. The three `type`
+# families produce three different things from an identically shaped result: a
+# single linear expression. Logistic regression gives a probability, SVM
+# classification gives a decision value whose sign picks the class, and SVR
+# gives a plain number. Cutting an SVM decision value at 0.5 as though it were
+# a probability misclassifies every row whose value lies between 0 and 0.5.
+
+#' @export
+tidypredict_output_type.LiblineaR <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (x$Type %in% liblinear_lr_types) {
+    return("prob")
+  }
+  if (x$Type %in% liblinear_svm_class_types) {
+    return("decision")
+  }
+  if (x$Type %in% liblinear_regression_types) {
+    return("numeric")
+  }
+
+  # parse_model_liblinear() rejects any other type, so reaching here means the
+  # two lists have drifted apart.
+  cli::cli_abort(
+    "Unsupported {.pkg LiblineaR} model {.arg type} {.val {x$Type}}.",
+    .internal = TRUE
+  )
+}
+
+#' @export
+tidypredict_outcome_levels.LiblineaR <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (x$Type %in% liblinear_regression_types) {
+    return(NULL)
+  }
+
+  # `ClassNames` is a factor whose element order reflects LiblineaR's internal
+  # ordering, not the outcome's level order. `parse_model_liblinear()` already
+  # normalises to the glm convention of predicting the second level, so report
+  # the levels rather than the elements to stay consistent with it.
+  levels(x$ClassNames)
+}
+
+#' @export
+tidypredict_normalized.LiblineaR <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # Binary only, and a single expression, so there is no set of per-level
+  # values to sum. `TRUE` would imply the caller can read probabilities for
+  # every level straight off the result, which it cannot.
+  NA
+}
+
 # Test --------------------------------------------
 
 #' @export
diff --git a/R/model-lightgbm.R b/R/model-lightgbm.R
index b96cefbd..0b0ac432 100644
--- a/R/model-lightgbm.R
+++ b/R/model-lightgbm.R
@@ -911,3 +911,42 @@ tidypredict_n_trees.lgb.Booster <- function(x, ...) {
   # trees actually returned rather than the number LightGBM reports.
   length(tidypredict_trees(x))
 }
+
+# Output metadata ---------------------------------
+
+# The same objective groups `build_fit_formula_lgb()` switches on: the
+# multiclass objectives softmax one raw score per class, the sigmoid objectives
+# give a single binary probability, and the rest stay numeric.
+lgb_parsed_objective <- function(x) {
+  x$general$params$objective %||% "regression"
+}
+
+#' @export
+tidypredict_output_type.pm_lgb <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  objective <- lgb_parsed_objective(x)
+  if (objective %in% c(lgb_multiclass_objectives, lgb_sigmoid_objectives)) {
+    return("prob")
+  }
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.pm_lgb <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # LightGBM is fit on integer labels. The multiclass expressions come back
+  # named `class_0`, `class_1` and so on, which are positions, not levels.
+  NULL
+}
+
+#' @export
+tidypredict_normalized.pm_lgb <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (lgb_parsed_objective(x) %in% lgb_multiclass_objectives) {
+    return(TRUE)
+  }
+  NA
+}
diff --git a/R/model-nnet.R b/R/model-nnet.R
index b457cb60..3d502108 100644
--- a/R/model-nnet.R
+++ b/R/model-nnet.R
@@ -224,3 +224,35 @@ tidypredict_test.nnet <- function(
     model$call
   )
 }
+
+# Output metadata ---------------------------------
+
+# `build_fit_formula_nnet()` returns one expression per level for a
+# classification fit, either from `expr_softmax()` or as `1 - p` and `p` for a
+# single logistic output unit, so the per-level values always sum to one. With
+# no levels recorded the single output unit is a plain numeric prediction.
+#' @export
+tidypredict_output_type.pm_nnet <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (is.null(parsed_model_classes(x))) {
+    return("numeric")
+  }
+  "prob"
+}
+
+#' @export
+tidypredict_outcome_levels.pm_nnet <- function(x, ...) {
+  rlang::check_dots_empty()
+  parsed_model_classes(x)
+}
+
+#' @export
+tidypredict_normalized.pm_nnet <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (is.null(parsed_model_classes(x))) {
+    return(NA)
+  }
+  TRUE
+}
diff --git a/R/model-partykit.R b/R/model-partykit.R
index 548f1697..abbc2263 100644
--- a/R/model-partykit.R
+++ b/R/model-partykit.R
@@ -333,3 +333,39 @@ tidypredict_class_exprs.party <- function(x, ...) {
 build_tree_formula.pm_tree_party <- function(model) {
   generate_nested_case_when_tree(model$tree_info, missing = "na")
 }
+
+# Output metadata ---------------------------------
+
+# `party_tree_info()` reads the mode off the response column of `fitted`, and
+# the parsed form does not record it, so the fitted object answers.
+party_is_classification <- function(x) {
+  is.factor(x$fitted[["(response)"]])
+}
+
+#' @export
+tidypredict_output_type.party <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (party_is_classification(x)) {
+    return("class")
+  }
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.party <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (party_is_classification(x)) {
+    return(levels(x$fitted[["(response)"]]))
+  }
+  NULL
+}
+
+#' @export
+tidypredict_normalized.party <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # A single expression, so there are no per-level values to sum.
+  NA
+}
diff --git a/R/model-rpart.R b/R/model-rpart.R
index 0297dc24..29a0a42f 100644
--- a/R/model-rpart.R
+++ b/R/model-rpart.R
@@ -287,3 +287,36 @@ rpart_classprob_tree_info <- function(model) {
 build_tree_formula.pm_tree_rpart <- function(model) {
   build_tree_formula_single(model)
 }
+
+# Output metadata ---------------------------------
+
+# The parsed form keeps only the tree, so the mode has to come off the fitted
+# object. `method` is the same signal `tidypredict_test.rpart()` uses.
+#' @export
+tidypredict_output_type.rpart <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (identical(x$method, "class")) {
+    return("class")
+  }
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.rpart <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (identical(x$method, "class")) {
+    return(attr(x, "ylevels"))
+  }
+  NULL
+}
+
+#' @export
+tidypredict_normalized.rpart <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # One expression, a class label for a classification tree and a number for a
+  # regression one, so there are no per-level values to sum either way.
+  NA
+}
diff --git a/R/model-rq.R b/R/model-rq.R
index 81f6c5bd..e912843c 100644
--- a/R/model-rq.R
+++ b/R/model-rq.R
@@ -50,3 +50,29 @@ split_rqs <- function(model) {
     }
   )
 }
+
+# Output metadata ---------------------------------
+
+# `parse_model.rqs()` returns a bare list of parsed models rather than a single
+# parsed model, so the default cannot route through it. Several `tau` give a
+# named list of expressions that looks exactly like a multiclass probability
+# list but holds unrelated quantile predictions.
+#' @export
+tidypredict_output_type.rqs <- function(x, ...) {
+  rlang::check_dots_empty()
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.rqs <- function(x, ...) {
+  rlang::check_dots_empty()
+  NULL
+}
+
+#' @export
+tidypredict_normalized.rqs <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # The names are quantiles, not levels, and the values do not sum to anything.
+  NA
+}
diff --git a/R/model-xgboost.R b/R/model-xgboost.R
index e2506dd3..9db18d59 100644
--- a/R/model-xgboost.R
+++ b/R/model-xgboost.R
@@ -601,3 +601,44 @@ tidypredict_n_trees.xgb.Booster <- function(x, ...) {
 
   length(tidypredict_trees(x))
 }
+
+# Output metadata ---------------------------------
+
+# `apply_xgb_objective()` decides this. The two logistic objectives get wrapped
+# in a logistic and so are probabilities. `binary:hinge` gets wrapped in
+# `as.numeric(score >= 0)`, which takes only the values 0 and 1: a hard class
+# prediction rather than a number, even though it is numerically typed.
+# Everything else stays a raw score on the response scale, and the multiclass
+# objectives are rejected outright.
+#' @export
+tidypredict_output_type.pm_xgb <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  objective <- x$general$params$objective
+  if (
+    identical(objective, "binary:logistic") ||
+      identical(objective, "reg:logistic")
+  ) {
+    return("prob")
+  }
+  if (identical(objective, "binary:hinge")) {
+    return("class")
+  }
+  "numeric"
+}
+
+#' @export
+tidypredict_outcome_levels.pm_xgb <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # xgboost is fit on a numeric label, so no fit ever records outcome levels.
+  NULL
+}
+
+#' @export
+tidypredict_normalized.pm_xgb <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # Multiclass objectives are unsupported, so the fit is always one expression.
+  NA
+}
diff --git a/R/output-metadata.R b/R/output-metadata.R
new file mode 100644
index 00000000..2ff31f77
--- /dev/null
+++ b/R/output-metadata.R
@@ -0,0 +1,343 @@
+#' Describe what a model's fitted expressions compute
+#'
+#' @description
+#' `tidypredict_fit()` returns expressions, but not what those expressions
+#' mean. A single expression could be a numeric prediction, a probability, or
+#' an uncalibrated decision value, and the three call for different handling
+#' downstream. These generics answer that question, so that a package
+#' generating code from the result does not have to keep its own list of which
+#' backend produces which shape.
+#'
+#' The metadata is asked of the *model*, not of the fitted expressions,
+#' deliberately. Attributes on the result do not survive the subsetting,
+#' `lapply()` and `unlist()` that callers apply to a multiclass result, which
+#' is the case that most needs describing.
+#'
+#' @param x A fitted model object.
+#' @param ... Additional arguments passed to methods.
+#'
+#' @returns
+#' `tidypredict_output_type()` returns a single string, one of:
+#'
+#' \describe{
+#'   \item{`"numeric"`}{A numeric prediction. `tidypredict_fit()` returns one
+#'     expression, or a named list of them for a multivariate outcome or for a
+#'     quantile regression with several `tau`.}
+#'   \item{`"prob"`}{A probability. Either one expression giving the
+#'     probability of the second outcome level, for a binary model, or a list
+#'     with one expression per level.}
+#'   \item{`"decision"`}{An uncalibrated decision value whose *sign* selects
+#'     the class. Not a probability, and not comparable to one: the cut is at
+#'     0, not 0.5.}
+#'   \item{`"class"`}{A hard class prediction, with no probability available.
+#'     Usually the class label as a string, but `xgboost`'s `binary:hinge`
+#'     objective gives a 0/1 indicator instead. What makes it `"class"` rather
+#'     than `"numeric"` is that only the class values can occur, so using it as
+#'     a numeric prediction is a mistake even when its type is numeric.}
+#' }
+#'
+#' `tidypredict_outcome_levels()` returns a character vector of outcome levels
+#' in model order, or `NULL`.
+#'
+#' `NULL` means two different things, and `tidypredict_output_type()`
+#' distinguishes them. For a `"numeric"` model it means there are no levels.
+#' For a `"prob"` or `"class"` model it means the fitted model **did not retain
+#' the outcome levels**, so any names on the result are positional
+#' placeholders and the caller has to supply the real levels from elsewhere.
+#' LightGBM and CatBoost multiclass models are in this position: they store
+#' integer labels and their expressions come back named `class_0`, `class_1`
+#' and so on.
+#'
+#' `tidypredict_normalized()` returns `TRUE` if the per-level values already
+#' sum to one across levels, `FALSE` if the caller has to normalize them, and
+#' `NA` when there are no per-level values to sum, which includes every
+#' single-expression model.
+#'
+#' At present no backend returns `FALSE`: every multiclass probability list
+#' goes through one shared softmax, so the values are always normalized
+#' already. The generic exists so that a caller can rely on that rather than
+#' having to know it, and so a future backend that does not normalize can say
+#' so instead of silently breaking the assumption.
+#'
+#' @details
+#' None of this is recoverable from the shape of the result, which is the whole
+#' reason for recording it. Two concrete cases:
+#'
+#' A binary `"prob"` model and a `"decision"` model both return exactly one
+#' expression. `LiblineaR` produces either, depending only on its `type`
+#' argument. Treating a decision value as a probability and cutting it at 0.5
+#' gives silently wrong classes for every row whose value falls between 0 and
+#' 0.5.
+#'
+#' A multiclass `"prob"` model and a `quantreg::rq()` fit with several `tau`
+#' both return a named list of expressions of the same length and structure.
+#' In the first the values sum to one across the list; in the second they are
+#' unrelated numeric predictions.
+#'
+#' @examplesIf rlang::is_installed("MASS")
+#' model <- lm(mpg ~ wt, data = mtcars)
+#' tidypredict_output_type(model)
+#' tidypredict_outcome_levels(model)
+#' tidypredict_normalized(model)
+#'
+#' @name tidypredict_metadata
+NULL
+
+#' @rdname tidypredict_metadata
+#' @export
+tidypredict_output_type <- function(x, ...) {
+  UseMethod("tidypredict_output_type")
+}
+
+# Most models describe themselves adequately once parsed, and the parsed
+# classes are far fewer than the model classes, so the default routes through
+# `parse_model()` the same way `tidypredict_fit.default()` does. A model whose
+# parsed form is not enough, such as LiblineaR, gets its own method.
+#' @export
+tidypredict_output_type.default <- function(x, ...) {
+  metadata_via_parsed(
+    x,
+    tidypredict_output_type,
+    "tidypredict_output_type",
+    ...
+  )
+}
+
+metadata_via_parsed <- function(
+  x,
+  generic,
+  generic_name,
+  ...,
+  call = rlang::caller_env()
+) {
+  # A parsed model reaching the default means its own type has no method;
+  # parsing again would recurse forever.
+  if (inherits(x, "parsed_model")) {
+    abort_no_metadata(x, generic_name, call = call)
+  }
+
+  has_parser <- any(map_lgl(
+    class(x),
+    ~ !is.null(utils::getS3method("parse_model", .x, optional = TRUE))
+  ))
+  if (!has_parser) {
+    abort_no_metadata(x, generic_name, call = call)
+  }
+
+  generic(parse_model(x), ...)
+}
+
+#' @rdname tidypredict_metadata
+#' @export
+tidypredict_outcome_levels <- function(x, ...) {
+  UseMethod("tidypredict_outcome_levels")
+}
+
+#' @export
+tidypredict_outcome_levels.default <- function(x, ...) {
+  metadata_via_parsed(
+    x,
+    tidypredict_outcome_levels,
+    "tidypredict_outcome_levels",
+    ...
+  )
+}
+
+#' @rdname tidypredict_metadata
+#' @export
+tidypredict_normalized <- function(x, ...) {
+  UseMethod("tidypredict_normalized")
+}
+
+#' @export
+tidypredict_normalized.default <- function(x, ...) {
+  metadata_via_parsed(x, tidypredict_normalized, "tidypredict_normalized", ...)
+}
+
+# Parsed model methods ----------------------------
+
+# `build_fit_formula()` applies the inverse link when `is_glm` is 1, so a
+# binomial model's expression already is a probability. Every other family
+# returns a conditional mean on the response scale, which is a number rather
+# than a probability: a Poisson count of 3.2 is not 320% of anything.
+#' @export
+tidypredict_output_type.pm_regression <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  if (
+    identical(x$general$is_glm, 1) && identical(x$general$family, "binomial")
+  ) {
+    return("prob")
+  }
+
+  "numeric"
+}
+
+#' @export
+tidypredict_normalized.pm_regression <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  # One expression, so there is no set of per-level values to sum.
+  NA
+}
+
+# The parsed form of a regression keeps coefficients, not the outcome. The
+# binary probability backends that do record their response levels (glm,
+# `lognet`, `ksvm`) answer on the fitted class instead.
+#' @export
+tidypredict_outcome_levels.pm_regression <- function(x, ...) {
+  rlang::check_dots_empty()
+  NULL
+}
+
+# Every one of these produces one expression per outcome level.
+#' @export
+tidypredict_output_type.pm_multiclass_regression <- function(x, ...) {
+  rlang::check_dots_empty()
+  "prob"
+}
+
+#' @export
+tidypredict_output_type.pm_naive_bayes <- function(x, ...) {
+  rlang::check_dots_empty()
+  "prob"
+}
+
+#' @export
+tidypredict_output_type.pm_nullmodel_classification <- function(x, ...) {
+  rlang::check_dots_empty()
+  "prob"
+}
+
+# All three keep the outcome levels in `classes` and name the expressions after
+# them, and all three finish with `expr_softmax()` or an explicit set of
+# probabilities, so the per-level values already sum to one.
+#' @export
+tidypredict_outcome_levels.pm_multiclass_regression <- function(x, ...) {
+  rlang::check_dots_empty()
+  parsed_model_classes(x)
+}
+
+#' @export
+tidypredict_normalized.pm_multiclass_regression <- function(x, ...) {
+  rlang::check_dots_empty()
+  TRUE
+}
+
+#' @export
+tidypredict_outcome_levels.pm_naive_bayes <- function(x, ...) {
+  rlang::check_dots_empty()
+  parsed_model_classes(x)
+}
+
+#' @export
+tidypredict_normalized.pm_naive_bayes <- function(x, ...) {
+  rlang::check_dots_empty()
+  TRUE
+}
+
+#' @export
+tidypredict_outcome_levels.pm_nullmodel_classification <- function(x, ...) {
+  rlang::check_dots_empty()
+  parsed_model_classes(x)
+}
+
+#' @export
+tidypredict_normalized.pm_nullmodel_classification <- function(x, ...) {
+  rlang::check_dots_empty()
+  TRUE
+}
+
+# `classes` is a bare character vector for some backends and a list of strings
+# for others, which is also how `build_fit_formula_*()` reads it.
+parsed_model_classes <- function(x) {
+  classes <- as.character(unlist(x$classes))
+  if (length(classes) == 0) {
+    return(NULL)
+  }
+  classes
+}
+
+# `pm_tree` spans ten backends with three different modes between them, and
+# `general$model` is what the rest of the tree code switches on, so switch on it
+# here too. `ranger`, `randomForest`, `cforest`, `aorsf`, `blackboost` and
+# `cubist` all abort at parse time for a classification fit, so a parsed model
+# under those names is necessarily a regression. `bagger` records its levels in
+# `general$classes` when, and only when, it is a classification.
+#
+# `rpart` and `party` are the exception: their parsed form keeps the tree and
+# nothing that says which mode it came from, so they answer on the fitted class
+# and a parsed model on its own cannot be described.
+tree_regression_models <- c(
+  "ranger",
+  "randomForest",
+  "cforest",
+  "aorsf",
+  "blackboost",
+  "cubist"
+)
+
+#' @export
+tidypredict_output_type.pm_tree <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  model <- x$general$model
+  if (model %in% tree_regression_models) {
+    return("numeric")
+  }
+  if (identical(model, "C5.0")) {
+    return("class")
+  }
+  if (identical(model, "bagger")) {
+    if (is.null(parsed_model_classes(x$general))) {
+      return("numeric")
+    }
+    return("class")
+  }
+
+  abort_no_metadata(x, "tidypredict_output_type")
+}
+
+#' @export
+tidypredict_outcome_levels.pm_tree <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  model <- x$general$model
+  # `parse_model.C5.0()` keeps the tree but not `model$levels`, so the parsed
+  # form genuinely has no levels to report. `tidypredict_outcome_levels.C5.0()`
+  # answers from the fitted object.
+  if (model %in% c(tree_regression_models, "C5.0")) {
+    return(NULL)
+  }
+  if (identical(model, "bagger")) {
+    return(parsed_model_classes(x$general))
+  }
+
+  abort_no_metadata(x, "tidypredict_outcome_levels")
+}
+
+#' @export
+tidypredict_normalized.pm_tree <- function(x, ...) {
+  rlang::check_dots_empty()
+
+  model <- x$general$model
+  if (model %in% c(tree_regression_models, "C5.0", "bagger")) {
+    # Every tree backend produces a single expression, a number or a class
+    # label, so there are never per-level values to sum.
+    return(NA)
+  }
+
+  abort_no_metadata(x, "tidypredict_normalized")
+}
+
+# Aborting rather than guessing "numeric". A wrong answer here is not a visible
+# failure downstream, it is mis-named or mis-thresholded predictions, so a
+# backend with no method must say so rather than be assumed benign.
+abort_no_metadata <- function(x, generic, call = rlang::caller_env()) {
+  cli::cli_abort(
+    "{.fn {generic}} is not available for models of class
+     {.cls {class(x)[[1]]}}.",
+    class = "tidypredict_no_metadata",
+    call = call
+  )
+}
diff --git a/man/tidypredict_metadata.Rd b/man/tidypredict_metadata.Rd
new file mode 100644
index 00000000..622af453
--- /dev/null
+++ b/man/tidypredict_metadata.Rd
@@ -0,0 +1,99 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/output-metadata.R
+\name{tidypredict_metadata}
+\alias{tidypredict_metadata}
+\alias{tidypredict_output_type}
+\alias{tidypredict_outcome_levels}
+\alias{tidypredict_normalized}
+\title{Describe what a model's fitted expressions compute}
+\usage{
+tidypredict_output_type(x, ...)
+
+tidypredict_outcome_levels(x, ...)
+
+tidypredict_normalized(x, ...)
+}
+\arguments{
+\item{x}{A fitted model object.}
+
+\item{...}{Additional arguments passed to methods.}
+}
+\value{
+\code{tidypredict_output_type()} returns a single string, one of:
+
+\describe{
+\item{\code{"numeric"}}{A numeric prediction. \code{tidypredict_fit()} returns one
+expression, or a named list of them for a multivariate outcome or for a
+quantile regression with several \code{tau}.}
+\item{\code{"prob"}}{A probability. Either one expression giving the
+probability of the second outcome level, for a binary model, or a list
+with one expression per level.}
+\item{\code{"decision"}}{An uncalibrated decision value whose \emph{sign} selects
+the class. Not a probability, and not comparable to one: the cut is at
+0, not 0.5.}
+\item{\code{"class"}}{A hard class prediction, with no probability available.
+Usually the class label as a string, but \code{xgboost}'s \code{binary:hinge}
+objective gives a 0/1 indicator instead. What makes it \code{"class"} rather
+than \code{"numeric"} is that only the class values can occur, so using it as
+a numeric prediction is a mistake even when its type is numeric.}
+}
+
+\code{tidypredict_outcome_levels()} returns a character vector of outcome levels
+in model order, or \code{NULL}.
+
+\code{NULL} means two different things, and \code{tidypredict_output_type()}
+distinguishes them. For a \code{"numeric"} model it means there are no levels.
+For a \code{"prob"} or \code{"class"} model it means the fitted model \strong{did not retain
+the outcome levels}, so any names on the result are positional
+placeholders and the caller has to supply the real levels from elsewhere.
+LightGBM and CatBoost multiclass models are in this position: they store
+integer labels and their expressions come back named \code{class_0}, \code{class_1}
+and so on.
+
+\code{tidypredict_normalized()} returns \code{TRUE} if the per-level values already
+sum to one across levels, \code{FALSE} if the caller has to normalize them, and
+\code{NA} when there are no per-level values to sum, which includes every
+single-expression model.
+
+At present no backend returns \code{FALSE}: every multiclass probability list
+goes through one shared softmax, so the values are always normalized
+already. The generic exists so that a caller can rely on that rather than
+having to know it, and so a future backend that does not normalize can say
+so instead of silently breaking the assumption.
+}
+\description{
+\code{tidypredict_fit()} returns expressions, but not what those expressions
+mean. A single expression could be a numeric prediction, a probability, or
+an uncalibrated decision value, and the three call for different handling
+downstream. These generics answer that question, so that a package
+generating code from the result does not have to keep its own list of which
+backend produces which shape.
+
+The metadata is asked of the \emph{model}, not of the fitted expressions,
+deliberately. Attributes on the result do not survive the subsetting,
+\code{lapply()} and \code{unlist()} that callers apply to a multiclass result, which
+is the case that most needs describing.
+}
+\details{
+None of this is recoverable from the shape of the result, which is the whole
+reason for recording it. Two concrete cases:
+
+A binary \code{"prob"} model and a \code{"decision"} model both return exactly one
+expression. \code{LiblineaR} produces either, depending only on its \code{type}
+argument. Treating a decision value as a probability and cutting it at 0.5
+gives silently wrong classes for every row whose value falls between 0 and
+0.5.
+
+A multiclass \code{"prob"} model and a \code{quantreg::rq()} fit with several \code{tau}
+both return a named list of expressions of the same length and structure.
+In the first the values sum to one across the list; in the second they are
+unrelated numeric predictions.
+}
+\examples{
+\dontshow{if (rlang::is_installed("MASS")) withAutoprint(\{ # examplesIf}
+model <- lm(mpg ~ wt, data = mtcars)
+tidypredict_output_type(model)
+tidypredict_outcome_levels(model)
+tidypredict_normalized(model)
+\dontshow{\}) # examplesIf}
+}
diff --git a/tests/testthat/_snaps/output-metadata.md b/tests/testthat/_snaps/output-metadata.md
new file mode 100644
index 00000000..c09f7c0c
--- /dev/null
+++ b/tests/testthat/_snaps/output-metadata.md
@@ -0,0 +1,16 @@
+# models with no metadata method say so
+
+    Code
+      tidypredict_output_type(structure(list(), class = "made_up_model"))
+    Condition
+      Error in `tidypredict_output_type()`:
+      ! `tidypredict_output_type()` is not available for models of class <made_up_model>.
+
+# an rpart parsed model cannot say which mode it came from
+
+    Code
+      tidypredict_output_type(pm)
+    Condition
+      Error in `tidypredict_output_type()`:
+      ! `tidypredict_output_type()` is not available for models of class <parsed_model>.
+
diff --git a/tests/testthat/test-output-metadata.R b/tests/testthat/test-output-metadata.R
new file mode 100644
index 00000000..f1458d8b
--- /dev/null
+++ b/tests/testthat/test-output-metadata.R
@@ -0,0 +1,577 @@
+test_that("regression models report numeric output", {
+  expect_identical(tidypredict_output_type(lm(mpg ~ wt, mtcars)), "numeric")
+  expect_identical(
+    tidypredict_output_type(glm(mpg ~ wt, data = mtcars)),
+    "numeric"
+  )
+})
+
+test_that("only binomial glms report prob, not every glm", {
+  binom <- glm(vs ~ wt, data = mtcars, family = binomial)
+  pois <- glm(carb ~ wt, data = mtcars, family = poisson)
+
+  expect_identical(tidypredict_output_type(binom), "prob")
+  expect_identical(tidypredict_output_type(pois), "numeric")
+
+  # The claim is about values, not labels: a binomial fit lands in [0, 1] and
+  # a Poisson one does not, so calling both "prob" would be wrong.
+  binom_vals <- rlang::eval_tidy(tidypredict_fit(binom), mtcars)
+  pois_vals <- rlang::eval_tidy(tidypredict_fit(pois), mtcars)
+
+  expect_true(all(binom_vals >= 0 & binom_vals <= 1))
+  expect_false(all(pois_vals >= 0 & pois_vals <= 1))
+})
+
+test_that("a single expression has no per-level sum to normalize", {
+  expect_identical(tidypredict_normalized(lm(mpg ~ wt, mtcars)), NA)
+})
+
+test_that("LiblineaR distinguishes probability from decision value", {
+  skip_if_not_installed("LiblineaR")
+
+  df <- iris[iris$Species != "virginica", ]
+  df$Species <- droplevels(df$Species)
+  x <- as.matrix(df[, 1:4])
+
+  lr <- LiblineaR::LiblineaR(data = x, target = df$Species, type = 0)
+  svm <- LiblineaR::LiblineaR(data = x, target = df$Species, type = 1)
+
+  expect_identical(tidypredict_output_type(lr), "prob")
+  expect_identical(tidypredict_output_type(svm), "decision")
+
+  # Same shape, one expression each, but different meaning. The logistic fit
+  # lies in [0, 1]; the SVM decision value straddles 0, so thresholding it at
+  # 0.5 would misclassify every row between 0 and 0.5.
+  lr_vals <- rlang::eval_tidy(tidypredict_fit(lr), df)
+  svm_vals <- rlang::eval_tidy(tidypredict_fit(svm), df)
+
+  expect_true(all(lr_vals >= 0 & lr_vals <= 1))
+  expect_false(all(svm_vals >= 0 & svm_vals <= 1))
+  expect_true(any(svm_vals < 0))
+})
+
+test_that("LiblineaR reports outcome levels for classifiers only", {
+  skip_if_not_installed("LiblineaR")
+
+  df <- iris[iris$Species != "virginica", ]
+  df$Species <- droplevels(df$Species)
+  x <- as.matrix(df[, 1:4])
+
+  svm <- LiblineaR::LiblineaR(data = x, target = df$Species, type = 1)
+  expect_identical(
+    tidypredict_outcome_levels(svm),
+    c("setosa", "versicolor")
+  )
+
+  # svr_eps is passed explicitly only to silence LiblineaR's default notice.
+  svr <- LiblineaR::LiblineaR(
+    data = x,
+    target = df$Sepal.Length,
+    type = 11,
+    svr_eps = 0.1
+  )
+  expect_identical(tidypredict_output_type(svr), "numeric")
+  expect_null(tidypredict_outcome_levels(svr))
+})
+
+test_that("models with no metadata method say so", {
+  expect_snapshot(
+    error = TRUE,
+    tidypredict_output_type(structure(list(), class = "made_up_model"))
+  )
+
+  expect_error(
+    tidypredict_output_type(structure(list(), class = "made_up_model")),
+    class = "tidypredict_no_metadata"
+  )
+})
+
+test_that("a parsed model with no method does not recurse", {
+  pm <- as_parsed_model(list(general = list(type = "made_up")))
+
+  expect_error(
+    tidypredict_output_type(pm),
+    class = "tidypredict_no_metadata"
+  )
+})
+
+# Value-level semantics ---------------------------
+
+eval_metadata_fit <- function(fit, df) {
+  if (!is.list(fit)) {
+    return(rlang::eval_tidy(fit, df))
+  }
+  do.call(cbind, lapply(fit, function(e) rlang::eval_tidy(e, df)))
+}
+
+test_that("a prob list really does sum to one across levels", {
+  skip_if_not_installed("nnet")
+
+  model <- nnet::multinom(Species ~ ., data = iris, trace = FALSE)
+
+  expect_identical(tidypredict_output_type(model), "prob")
+  expect_identical(tidypredict_normalized(model), TRUE)
+  expect_identical(tidypredict_outcome_levels(model), levels(iris$Species))
+
+  probs <- eval_metadata_fit(tidypredict_fit(model), iris)
+  expect_equal(rowSums(probs), rep(1, nrow(iris)), tolerance = 1e-8)
+})
+
+test_that("a numeric list of the same shape does not sum to one", {
+  skip_if_not_installed("quantreg")
+
+  model <- quantreg::rq(mpg ~ wt, tau = c(0.25, 0.5, 0.75), data = mtcars)
+
+  expect_identical(tidypredict_output_type(model), "numeric")
+  expect_null(tidypredict_outcome_levels(model))
+  expect_identical(tidypredict_normalized(model), NA)
+
+  # `parse_model_lm()` calls `summary()` on each single-quantile fit, which
+  # warns about a nonunique solution and is unrelated to the metadata.
+  fit <- suppressWarnings(tidypredict_fit(model))
+  quantiles <- eval_metadata_fit(fit, mtcars)
+  expect_length(colnames(quantiles), 3)
+  # The point of the metadata: same shape as a multiclass probability list,
+  # so treating it as one would be wrong.
+  expect_false(any(abs(rowSums(quantiles) - 1) < 1))
+})
+
+test_that("a class model returns labels drawn from its levels", {
+  skip_if_not_installed("rpart")
+
+  model <- rpart::rpart(Species ~ ., data = iris, method = "class")
+  levs <- tidypredict_outcome_levels(model)
+
+  expect_identical(tidypredict_output_type(model), "class")
+  expect_identical(levs, levels(iris$Species))
+  expect_identical(tidypredict_normalized(model), NA)
+
+  predicted <- eval_metadata_fit(tidypredict_fit(model), iris)
+  expect_type(predicted, "character")
+  expect_in(unique(predicted), levs)
+})
+
+test_that("a single prob expression stays inside [0, 1]", {
+  skip_if_not_installed("glmnet")
+
+  df <- iris[iris$Species != "virginica", ]
+  df$Species <- droplevels(df$Species)
+  model <- glmnet::glmnet(
+    as.matrix(df[, 1:4]),
+    df$Species,
+    family = "binomial",
+    lambda = 0.01
+  )
+
+  expect_identical(tidypredict_output_type(model), "prob")
+  expect_identical(tidypredict_outcome_levels(model), levels(df$Species))
+  expect_identical(tidypredict_normalized(model), NA)
+
+  probs <- eval_metadata_fit(tidypredict_fit(model), df)
+  expect_all_true(probs >= 0 & probs <= 1)
+})
+
+# Linear and additive backends --------------------
+
+test_that("glm only reports levels for a factor binomial outcome", {
+  df <- iris[iris$Species != "virginica", ]
+  df$Species <- droplevels(df$Species)
+
+  factor_fit <- glm(Species ~ Sepal.Length, data = df, family = binomial)
+  numeric_fit <- glm(vs ~ wt, data = mtcars, family = binomial)
+
+  expect_identical(tidypredict_outcome_levels(factor_fit), levels(df$Species))
+  # Fit on a 0/1 numeric, so the model kept no levels to report.
+  expect_null(tidypredict_outcome_levels(numeric_fit))
+  expect_null(tidypredict_outcome_levels(glm(mpg ~ wt, data = mtcars)))
+})
+
+test_that("glmnet families report the right type", {
+  skip_if_not_installed("glmnet")
+
+  gaussian <- glmnet::glmnet(mtcars[, -1], mtcars$mpg, lambda = 1)
+  multinomial <- glmnet::glmnet(
+    as.matrix(iris[, 1:4]),
+    iris$Species,
+    family = "multinomial",
+    lambda = 0.01
+  )
+
+  expect_identical(tidypredict_output_type(gaussian), "numeric")
+  expect_null(tidypredict_outcome_levels(gaussian))
+
+  expect_identical(tidypredict_output_type(multinomial), "prob")
+  expect_identical(
+    tidypredict_outcome_levels(multinomial),
+    levels(iris$Species)
+  )
+  expect_identical(tidypredict_normalized(multinomial), TRUE)
+})
+
+test_that("earth and xrf report numeric", {
+  skip_if_not_installed("earth")
+
+  expect_identical(
+    tidypredict_output_type(earth::earth(mpg ~ ., data = mtcars)),
+    "numeric"
+  )
+})
+
+test_that("ksvm distinguishes regression from binary classification", {
+  skip_if_not_installed("kernlab")
+
+  df <- iris[iris$Species != "virginica", ]
+  df$Species <- droplevels(df$Species)
+
+  svr <- kernlab::ksvm(
+    mpg ~ .,
+    data = mtcars,
+    kernel = "vanilladot",
+    type = "eps-svr"
+  )
+  svc <- kernlab::ksvm(
+    Species ~ .,
+    data = df,
+    kernel = "vanilladot",
+    prob.model = TRUE
+  )
+
+  expect_identical(tidypredict_output_type(svr), "numeric")
+  # `lev()` holds the sorted response values for an SVR, which are not levels.
+  expect_null(tidypredict_outcome_levels(svr))
+
+  expect_identical(tidypredict_output_type(svc), "prob")
+  expect_identical(tidypredict_outcome_levels(svc), levels(df$Species))
+  expect_identical(tidypredict_normalized(svc), NA)
+})
+
+test_that("nullmodel reports its mode", {
+  skip_if_not_installed("parsnip")
+
+  reg <- parsnip::nullmodel(mtcars[, -1], mtcars$mpg)
+  cls <- parsnip::nullmodel(iris[, 1:4], iris$Species)
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_null(tidypredict_outcome_levels(reg))
+  expect_identical(tidypredict_normalized(reg), NA)
+
+  expect_identical(tidypredict_output_type(cls), "prob")
+  expect_identical(tidypredict_outcome_levels(cls), levels(iris$Species))
+  expect_identical(tidypredict_normalized(cls), TRUE)
+})
+
+test_that("mixOmics separates regression from discriminant analysis", {
+  skip_if_not_installed("mixOmics")
+
+  pls <- mixOmics::pls(as.matrix(iris[, 1:3]), iris[, 4], ncomp = 2)
+  plsda <- mixOmics::plsda(as.matrix(iris[, 1:4]), iris$Species, ncomp = 2)
+
+  expect_identical(tidypredict_output_type(pls), "numeric")
+  expect_null(tidypredict_outcome_levels(pls))
+
+  expect_identical(tidypredict_output_type(plsda), "prob")
+  expect_identical(tidypredict_outcome_levels(plsda), levels(iris$Species))
+  expect_identical(tidypredict_normalized(plsda), TRUE)
+})
+
+# Discriminant and naive Bayes backends -----------
+
+test_that("discriminant backends report a normalized prob list", {
+  skip_if_not_installed("MASS")
+  skip_if_not_installed("mda")
+  skip_if_not_installed("sda")
+  skip_if_not_installed("sparsediscrim")
+
+  models <- list(
+    lda = MASS::lda(Species ~ ., data = iris),
+    qda = MASS::qda(Species ~ ., data = iris),
+    fda = mda::fda(Species ~ ., data = iris),
+    sda = sda::sda(as.matrix(iris[, 1:4]), iris$Species, verbose = FALSE),
+    lda_diag = sparsediscrim::lda_diag(Species ~ ., data = iris)
+  )
+
+  expect_all_equal(
+    vapply(models, tidypredict_output_type, character(1)),
+    "prob"
+  )
+  expect_all_true(vapply(models, tidypredict_normalized, logical(1)))
+  for (model in models) {
+    expect_identical(tidypredict_outcome_levels(model), levels(iris$Species))
+  }
+})
+
+test_that("naive Bayes backends report a normalized prob list", {
+  skip_if_not_installed("naivebayes")
+  skip_if_not_installed("klaR")
+
+  nb <- naivebayes::naive_bayes(Species ~ ., data = iris)
+  knb <- klaR::NaiveBayes(Species ~ ., data = iris)
+
+  expect_identical(tidypredict_output_type(nb), "prob")
+  expect_identical(tidypredict_output_type(knb), "prob")
+  expect_identical(tidypredict_outcome_levels(nb), levels(iris$Species))
+  expect_identical(tidypredict_outcome_levels(knb), levels(iris$Species))
+  expect_identical(tidypredict_normalized(nb), TRUE)
+  expect_identical(tidypredict_normalized(knb), TRUE)
+})
+
+test_that("nnet reports its mode", {
+  skip_if_not_installed("nnet")
+
+  cls <- nnet::nnet(Species ~ ., data = iris, size = 1, trace = FALSE)
+  reg <- nnet::nnet(
+    mpg ~ wt,
+    data = mtcars,
+    size = 1,
+    trace = FALSE,
+    linout = TRUE
+  )
+
+  expect_identical(tidypredict_output_type(cls), "prob")
+  expect_identical(tidypredict_outcome_levels(cls), levels(iris$Species))
+  expect_identical(tidypredict_normalized(cls), TRUE)
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_null(tidypredict_outcome_levels(reg))
+  expect_identical(tidypredict_normalized(reg), NA)
+})
+
+# Tree backends ----------------------------------
+
+test_that("regression-only forests report numeric", {
+  skip_if_not_installed("ranger")
+  skip_if_not_installed("randomForest")
+  skip_if_not_installed("Cubist")
+
+  models <- list(
+    ranger = ranger::ranger(mpg ~ ., data = mtcars, num.trees = 2),
+    randomForest = randomForest::randomForest(mpg ~ ., mtcars, ntree = 2),
+    cubist = Cubist::cubist(mtcars[, -1], mtcars$mpg)
+  )
+
+  expect_all_equal(
+    vapply(models, tidypredict_output_type, character(1)),
+    "numeric"
+  )
+  for (model in models) {
+    expect_null(tidypredict_outcome_levels(model))
+    expect_identical(tidypredict_normalized(model), NA)
+  }
+})
+
+test_that("rpart reports its mode from `method`", {
+  skip_if_not_installed("rpart")
+
+  reg <- rpart::rpart(mpg ~ wt, data = mtcars)
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_null(tidypredict_outcome_levels(reg))
+  expect_identical(tidypredict_normalized(reg), NA)
+})
+
+test_that("partykit ctree reports its mode from the response", {
+  skip_if_not_installed("partykit")
+
+  reg <- partykit::ctree(mpg ~ wt, data = mtcars)
+  cls <- partykit::ctree(Species ~ ., data = iris)
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_null(tidypredict_outcome_levels(reg))
+
+  expect_identical(tidypredict_output_type(cls), "class")
+  expect_identical(tidypredict_outcome_levels(cls), levels(iris$Species))
+  expect_identical(tidypredict_normalized(cls), NA)
+})
+
+test_that("C5.0 reports a class label", {
+  skip_if_not_installed("C50")
+
+  model <- C50::C5.0(iris[, 1:4], iris$Species, trials = 2)
+
+  expect_identical(tidypredict_output_type(model), "class")
+  expect_identical(tidypredict_outcome_levels(model), levels(iris$Species))
+  expect_identical(tidypredict_normalized(model), NA)
+})
+
+test_that("bagger reports its mode", {
+  skip_if_not_installed("baguette")
+
+  reg <- baguette::bagger(mpg ~ wt + cyl, data = mtcars, times = 2)
+  cls <- baguette::bagger(Species ~ ., data = iris, times = 2)
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_null(tidypredict_outcome_levels(reg))
+
+  expect_identical(tidypredict_output_type(cls), "class")
+  expect_identical(tidypredict_outcome_levels(cls), levels(iris$Species))
+  expect_identical(tidypredict_normalized(cls), NA)
+})
+
+test_that("aorsf and mboost report numeric", {
+  skip_if_not_installed("aorsf")
+  skip_if_not_installed("mboost")
+
+  orsf <- aorsf::orsf(mtcars, mpg ~ ., n_tree = 2)
+  bb <- mboost::blackboost(mpg ~ wt, data = mtcars)
+
+  expect_identical(tidypredict_output_type(orsf), "numeric")
+  expect_identical(tidypredict_output_type(bb), "numeric")
+  expect_null(tidypredict_outcome_levels(orsf))
+  expect_null(tidypredict_outcome_levels(bb))
+  expect_identical(tidypredict_normalized(orsf), NA)
+  expect_identical(tidypredict_normalized(bb), NA)
+})
+
+test_that("dbarts reports numeric", {
+  skip_if_not_installed("dbarts")
+
+  model <- dbarts::bart2(
+    mpg ~ wt,
+    data = mtcars,
+    n.trees = 5,
+    n.samples = 5,
+    n.burn = 5,
+    verbose = FALSE,
+    keepTrees = TRUE
+  )
+
+  expect_identical(tidypredict_output_type(model), "numeric")
+  expect_null(tidypredict_outcome_levels(model))
+  expect_identical(tidypredict_normalized(model), NA)
+})
+
+# Boosted backends -------------------------------
+
+test_that("xgboost reads its mode from the objective", {
+  skip_if_not_installed("xgboost")
+
+  reg <- xgboost::xgb.train(
+    list(objective = "reg:squarederror", max_depth = 2),
+    xgboost::xgb.DMatrix(as.matrix(mtcars[, -1]), label = mtcars$mpg),
+    nrounds = 2,
+    verbose = 0
+  )
+  binary <- xgboost::xgb.train(
+    list(objective = "binary:logistic", max_depth = 2),
+    xgboost::xgb.DMatrix(as.matrix(mtcars[, -8]), label = mtcars$vs),
+    nrounds = 2,
+    verbose = 0
+  )
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_identical(tidypredict_output_type(binary), "prob")
+  # xgboost is fit on a numeric label, so no fit records outcome levels.
+  expect_null(tidypredict_outcome_levels(binary))
+  expect_identical(tidypredict_normalized(binary), NA)
+})
+
+test_that("xgboost binary:hinge is a class prediction, not a number", {
+  skip_if_not_installed("xgboost")
+
+  hinge <- xgboost::xgb.train(
+    list(objective = "binary:hinge", max_depth = 2),
+    xgboost::xgb.DMatrix(as.matrix(mtcars[, -8]), label = mtcars$vs),
+    nrounds = 3,
+    verbose = 0
+  )
+
+  expect_identical(tidypredict_output_type(hinge), "class")
+
+  # The label follows from the values: `as.numeric(score >= 0)` can only ever
+  # be 0 or 1, so treating it as a numeric prediction would be wrong even
+  # though its type is numeric.
+  values <- rlang::eval_tidy(tidypredict_fit(hinge), mtcars)
+  expect_setequal(unique(values), c(0, 1))
+})
+
+test_that("lightgbm reads its mode from the objective", {
+  skip_if_not_installed("lightgbm")
+
+  reg <- lightgbm::lgb.train(
+    list(objective = "regression", num_leaves = 3, verbose = -1),
+    lightgbm::lgb.Dataset(as.matrix(mtcars[, -1]), label = mtcars$mpg),
+    nrounds = 2,
+    verbose = -1
+  )
+  binary <- lightgbm::lgb.train(
+    list(objective = "binary", num_leaves = 3, verbose = -1),
+    lightgbm::lgb.Dataset(as.matrix(mtcars[, -8]), label = mtcars$vs),
+    nrounds = 2,
+    verbose = -1
+  )
+  multi <- lightgbm::lgb.train(
+    list(objective = "multiclass", num_class = 3, num_leaves = 3, verbose = -1),
+    lightgbm::lgb.Dataset(
+      as.matrix(iris[, 1:4]),
+      label = as.integer(iris$Species) - 1
+    ),
+    nrounds = 2,
+    verbose = -1
+  )
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_identical(tidypredict_normalized(reg), NA)
+
+  expect_identical(tidypredict_output_type(binary), "prob")
+  expect_identical(tidypredict_normalized(binary), NA)
+
+  expect_identical(tidypredict_output_type(multi), "prob")
+  expect_identical(tidypredict_normalized(multi), TRUE)
+  # LightGBM stores integer labels, so the `class_*` names are positional and
+  # the caller has to supply the real levels.
+  expect_null(tidypredict_outcome_levels(multi))
+  expect_named(tidypredict_fit(multi), c("class_0", "class_1", "class_2"))
+})
+
+test_that("catboost reads its mode from the objective", {
+  skip_if_not_installed("catboost")
+
+  train <- function(X, y, loss_function) {
+    pool <- catboost_catboost.load_pool(
+      X,
+      label = y,
+      feature_names = as.list(colnames(X))
+    )
+    catboost_catboost.train(
+      pool,
+      params = list(
+        iterations = 5L,
+        depth = 2L,
+        learning_rate = 0.5,
+        loss_function = loss_function,
+        logging_level = "Silent",
+        allow_writing_files = FALSE
+      )
+    )
+  }
+
+  reg <- train(
+    data.matrix(mtcars[, c("cyl", "disp")]),
+    mtcars$mpg,
+    "RMSE"
+  )
+  multi <- train(
+    data.matrix(iris[, 1:4]),
+    as.integer(iris$Species) - 1L,
+    "MultiClass"
+  )
+
+  expect_identical(tidypredict_output_type(reg), "numeric")
+  expect_identical(tidypredict_normalized(reg), NA)
+  expect_null(tidypredict_outcome_levels(reg))
+
+  expect_identical(tidypredict_output_type(multi), "prob")
+  expect_identical(tidypredict_normalized(multi), TRUE)
+  # CatBoost stores integer labels, so the `class_*` names are positional.
+  expect_null(tidypredict_outcome_levels(multi))
+})
+
+# Parsed models ----------------------------------
+
+test_that("an rpart parsed model cannot say which mode it came from", {
+  skip_if_not_installed("rpart")
+
+  pm <- parse_model(rpart::rpart(Species ~ ., data = iris, method = "class"))
+
+  expect_snapshot(error = TRUE, tidypredict_output_type(pm))
+})