predict.ref_fit

Usage

predict.ref_fit(
  object,
  newdata,
  type = c("scores", "distribution", "harmonised"),
  uncertainty = c("total", "conditional"),
  outcomes = NULL,
  allow_extrapolation = FALSE,
  n_draw = NULL,
  by = NULL,
  to = NULL,
  ...
)

One method serves plain, adapted, calibrated, and frozen fits. The steps are: engine prediction, then site adaptation (shifting location and scale, including any coefficient draws), then the response transform of ref_spec() (which puts the predictive back on the scale of y), then scoring, then PIT recalibration, and finally extrapolation masking.

Arguments

Argument	Description
object	A ref_fit.
newdata	Data frame of target observations.
type	"scores", "distribution", or "harmonised". Distributions carry adaptation but not the calibration map, which acts on probabilities. "harmonised" returns newdata with each outcome replaced by its conditional quantile mapping out of the group effect: $y_h = Q(F(y \mid x, g) \mid x, g_0)$, where $F$ is the predictive for the row’s own group (after adaptation) and the target $g_0$ is either the model with the smooth terms of by excluded (a zero group effect; a parametric by term is set to its reference level) and no adaptation, or the level to. The map is deterministic and monotone, so harmonised values scored under the target predictive have exactly the original centiles; the other columns of newdata (including by) are returned unchanged.
uncertainty	"total" (default) integrates the scores over coefficient draws (or uses the analytic Gaussian total when the location is identity-linked with constant scale); "conditional" uses the point estimates.
outcomes	Optional subset of outcomes.
allow_extrapolation	If FALSE (default), rows with support "out" or "new_group" have NA for z, centile, tail_prob, tail_surprisal, and log_density (median and residual are kept). An unseen level of a grouping covariate is extrapolation unless the fit has been adapted to it with ref_adapt(); a parametric factor has no prediction at all for such a row.
n_draw	Number of coefficient draws for "total"; defaults to spec$control$n_draw or 200. Draws are seeded from spec$control$seed (default 1) so repeated calls agree exactly. With 200 draws the Monte Carlo error of a score is about 0.02 in z near the centre and larger in the far tails; it falls as $1/\sqrt{n_{draw}}$, so pass n_draw = 2000 (or set control = list(n_draw = 2000) in ref_spec()) for reporting extreme centiles. The constant-scale Gaussian uses the exact analytic total and is unaffected.
by	Grouping covariate for "harmonised"; defaults to the adaptation’s by or, failing that, the single factor entering a smooth (such as s(site, bs = "re")).
to	Optional target level for "harmonised"; NULL (default) removes the group effect.
...	Unused.

Value

For "scores", a ref_scores tibble with one row per observation and outcome. status is "ok", "missing_predictor" (a covariate is NA), "new_group" (an unseen factor level), or the fit status of a failed outcome. For "distribution", a tibble with .id and one distribution column per outcome (missing distributions for failed outcomes). For "harmonised", newdata with harmonised outcome columns.

Examples

ref <- ref_simulate(300, site_shift = c(0, 1, -1, 0), seed = 1)
spec <- ref_spec(ref_gaussian(), ~ s(age, k = 5) + sex + s(site, bs = "re"))
fit <- ref_fit(spec, ref, outcomes = "y")
# synthetic subjects from the reference: one draw per covariate-grid row
grid <- expand.grid(age = c(30, 50, 70), sex = factor("F"), site = factor("A"))
dists <- predict(fit, grid, type = "distribution")
synthetic <- data.frame(grid, y = unlist(distributional::generate(dists$y, 1)))
# site effects mapped out of the observed outcomes
harm <- predict(fit, ref, type = "harmonised")
tapply(harm$y - ref$y, ref$site, mean)

         A          B          C          D
 0.1971966 -0.9628965  0.8824625 -0.1164318