Fitted values

brms's fitted(): a summary of the expected response, in the columns Estimate, Est.Error, Q2.5 and Q97.5. The estimate is the modelled response at the estimates, conditional on the random-effect modes, which is frm_linpred(object, type = "response") on the training data for every family.

Usage

# S3 method for class 'frmtmb_fit'
fitted(
  object,
  newdata = NULL,
  re_formula = NULL,
  scale = c("response", "linear"),
  resp = NULL,
  dpar = NULL,
  nlpar = NULL,
  ndraws = NULL,
  draw_ids = NULL,
  sort = FALSE,
  summary = TRUE,
  robust = FALSE,
  probs = c(0.025, 0.975),
  ...,
  allow_new_levels = FALSE,
  sample_new_levels = NULL
)

Arguments

object: A frmtmb_fit.
newdata: Optional data frame to evaluate on. Defaults to the training data.
re_formula: Which group-level terms enter the answer: NULL (default) keeps all of them, NA keeps none, and a one-sided formula keeps the terms it names, as in brms (see frm_linpred()). brms's spelling, and the only one: lme4's re.form is not accepted here.
scale: "response" (default) for the modelled response, or "linear" for the linear predictor. brms's spelling of what predict.frmtmb_fit() calls type.
resp: For multivariate fits: the response or responses to report. NULL (default) reports all of them, as brms does.
dpar: Which distributional parameter to report instead of the mean.
nlpar: brms's name for a non-linear parameter, which is a distributional parameter here: a synonym for dpar, and giving both is an error.
ndraws: brms's number of posterior draws to summarize. It sets the Monte Carlo error of brms's answer, not the quantity answered, and a maximum-likelihood fit answers that quantity at its estimate with no Monte Carlo error. So it is accepted, changes nothing, and a message says so.
draw_ids, sort, summary, robust: brms's arguments, in brms's positions so that a positional brms call asks the same question. Each needs posterior draws, and a maximum-likelihood fit has none, so each is refused by name with the reason and with the place it does work: frmtmb.sample's posterior_epred(). The default of each is accepted and changes nothing.
probs: Probabilities of the two quantile columns. brms takes the posterior quantiles there; here they are the ends of the Wald interval at those probabilities.
...: Refused. An argument this method does not have is an error naming it, because a swallowed re_formula returned the conditional fit and said nothing.
allow_new_levels: Predict unseen grouping-factor levels at the population level instead of erroring, as in frm_linpred().
sample_new_levels: brms's argument. "gaussian", or NULL for it: an unseen level's effect is taken from its block's estimated covariance, so its variance widens Est.Error and the interval. "old_levels" gives each unseen level one seen level of its grouping factor, chosen at random once per call as brms chooses it (and as predict.frmtmb_fit() chooses it), and every term of that factor reads the same seen level: the row takes that level's fitted effects, and its Est.Error carries their conditional variance rather than the marginal one. One choice serves every response of a multivariate call. Two different unseen values of one mm() term read two seen levels; brms 2.23.0 gives them one, because it numbers unseen values per member rather than by value. It acts only with allow_new_levels = TRUE, as in brms. "uncertainty" chooses a level per posterior draw, and a maximum-likelihood fit has none, so it is refused by name.

Value

An n x 4 matrix with the columns Estimate, Est.Error and one per entry of probs. For an ordinal or categorical family an n x 4 x K array, the third dimension named P(Y = k), which is brms's shape; at scale = "linear", an ordinal fit with a cs() term gives n x 4 x (K - 1), one linear predictor per threshold, named eta1, eta2, ... as brms names them. For a multivariate fit asked for more than one response, an n x 4 x nresp array with the third dimension named by response, which is brms's shape too; each cell's Est.Error is the standard error of that response's own expected value. A category-valued response contributes one layer per category, named P(Y = k), in its place, as brms stacks it. The ROW dimnames are NULL, as brms's are; the data's row names are on frm_linpred() and model.frame().

Details

A maximum-likelihood fit has no draws to summarize, so Est.Error is the delta-method standard error frm_linpred() reports for the same quantity, and the Q columns are the Wald interval at those probabilities. The estimate is at the modes, and the standard error carries the uncertainty in them: at a grouping level the fit saw, the level's group effect enters with its covariance taken jointly with the fixed effects from the joint precision, the frequentist analogue of the posterior of that effect in brms's posterior_epred(). The interval covers the expected response at a known level with the nominal coverage averaged over the groups, not for one group's realized effect. The interval is around the EXPECTED response and carries no observation noise; predict.frmtmb_fit() is the predictive interval that does.

A nonlinear body's Est.Error is the delta method through the body's taped Jacobian. A predictor whose standard error this package cannot produce, such as a structured likelihood, reports NA in Est.Error and in the Q columns rather than a number it does not have. The estimate is unaffected.

The Q columns are the interval on the scale of the link, mapped through its inverse, so they stay in the support of the expected response, and they are not symmetric about the estimate unless the link is the identity.

New rows of a model with me()

A new row has no latent value, and brms's draws take one from N(x, sdx). The estimate here is the mean of the expected response over that draw, by Gauss-Hermite quadrature on the response scale. On a log link it is exp((b sdx)^2 / 2) times the value at x, and on an identity link it is that value. Est.Error holds the draw's variance beside the coefficients', so it is not coefficient uncertainty alone, and the Q columns are the draw's quantiles. frm_linpred() stays at the observed x. Two predictors that read one term read one latent value, a stated divergence from brms, and a category probability is a joint function of both, so its mean moves: 0.3336 against brms's 0.3106. A mixture and a mu with a zi that read one term do the same.

The product of two me() values on one row is read at first order: its mean is low by 7 percent at 0.8 and 0.6 on a log link.

A me() value that a nonlinear parameter reads is one of the draws below: fitted() integrates the body over it. frm_linpred(se.fit = TRUE) and frm_lp_basis() read it through the body at first order.

The draws of a new row

A new row can read three draws that are not coefficient uncertainty: a me() value, a latent residual, and with allow_new_levels = TRUE the effects of a grouping level the fit did not see. brms draws each of them, and its posterior_epred() averages over them; a new level is drawn as brms's sample_new_levels = "gaussian" draws it, not as its default "uncertainty". The estimate here is that mean. On a log link that reads the mean, and on a lognormal's mu, it is the closed form: exp(v / 2) times the value at the draws' mean of 0, with v their variance on the link scale, so at a new level of (1 | g) at sd 1 it is 1.65 times the median. Elsewhere it is Gauss-Hermite quadrature over the draws. The draws are independent of each other, and each one moves every predictor that reads it, so a latent on mu beside a me() value on hu is read in two dimensions; a draw that moves no cell of what is read costs nothing. sample_new_levels = "old_levels" reads a level the fit saw and takes no mean. scale = "linear" stays at the effects' mean, with their variance in Est.Error.

Est.Error is the draws' sd beside the coefficients'. At a new level it holds the block's estimated sd as known: the sd's own uncertainty is not in it, where brms's posterior carries it (0.90 of brms's on a poisson (1 | g) of 40 groups). With two or more directions the Q columns are the draw's quantiles along the direction in which the value moves most, widened to the whole variance, and covered 0.926 to 0.957 of the draws' own distribution on four hurdle rows. Past two directions the rule is checked against two fewer nodes, and a warning of class frmtmb_draw_rule gives the change when the sd moves by more than 2 percent. Past four directions that move the value, the four that move it most are kept, with a warning of class frmtmb_draw_dims.

One mean does not exist. A lognormal mean is exp(mu + sigma^2 / 2), and when sigma reads a draw on its log link, the exact mean over the draw is infinite at any size of the draw, from its far tail. This holds for a lognormal component's sigma<k> of a mixture and for a nonlinear parameter that a body for sigma reads. The estimate and Est.Error are then the 11-node Gauss-Hermite mean and sd, which reach 5.19 sds of the draw as a sampler's draws do, with one warning per call of class frmtmb_mean_not_finite that gives the largest change 9 nodes make. On a sigma ~ (1 | g) at sd 0.27 this reads 2.0897 against 2.0897 from 4e6 draws; at a draw sd of 0.4 on sigma, no mean of draws is stable, and neither is this; where a node overflows the row reads NaN. The Q columns are the draw's quantiles. Any other mean whose rule does not settle as it grows from 21 or 41 nodes is read the same way, with a warning of class frmtmb_mean_unsettled.

At a row whose mi() value left the objective, a reader's fit reads the value's mean. Under a nonlinear link this is below brms's average: a median of 25% for a log link on a Gamma value of shape 1.2.

A distributional parameter with no predictor

With scale = "linear", dpar reads a distributional parameter with no predictor on its link scale, as it reads one with a predictor. This diverges from brms 2.23.0, which returns the stored natural value. For one model, brms reads 7.48 with shape unmodeled and 2.01 with shape ~ 1, and frmtmb reads 2.01 in both. On another fit, frmtmb reads 0.776 where brms reads 2.17. scale = "response" reads the natural value in both packages.

Ordinal responses

The Q columns of a category probability are a Wald interval on its logit, as in conditional_effects(), so they stay inside 0 and 1. On a new me() row they are the draw's own quantiles, moved out on the logit by the coefficients' share of the variance.

An ordinal response has no mean, so fitted() summarizes the K category probabilities, with rows of Estimate summing to one, which is the brms fitted() convention. cs() terms are honored. The standard error of a category probability is the finite-difference delta method over the whole outer parameter vector, because the probability depends on the thresholds and the cs() coefficients as well as on the linear predictor. On a mixed ordinal fit the group effects join the differenced vector, with their covariance taken jointly with the parameters from the joint precision, so the standard error carries their uncertainty as the scalar route does. The estimates are unaffected. The latent linear predictor, which is where the coefficients live, is frm_linpred(object, type = "link").

Residual autocorrelation with `cov = FALSE`

Under ar(), ma() or arma() with cov = FALSE the mean of a row reads the residuals of the earlier rows of its group, taken from newdata's response. Where that response is NA, or absent, brms fills it with a draw from the family at the row's one-step mean, per posterior draw, so its fitted() carries the spread of the fill. A maximum likelihood fit has no draws to carry that spread through, so fitted() here fills the response with its expected value, the one-step mean itself. That is the mean of brms's fill at fixed parameters, for gaussian and student, the families cov = FALSE takes. Est.Error is then the parameter uncertainty alone, and is smaller than brms's on the rows after a missing response: this is a deliberate divergence. predict.frmtmb_fit() fills with draws, as brms does, and so does frmtmb.sample's posterior_epred().

Examples

set.seed(1)
dd <- data.frame(x = rnorm(100))
dd$y <- rpois(100, exp(0.3 + 0.4 * dd$x))
fit <- frm(bf(y ~ x) + poisson(), data = dd)
head(fitted(fit))
#>       Estimate Est.Error      Q2.5     Q97.5
#> [1,] 0.7751024 0.1130250 0.5824215 1.0315275
#> [2,] 1.1849807 0.1126043 0.9836142 1.4275712
#> [3,] 0.6946365 0.1132732 0.5046074 0.9562283
#> [4,] 2.4828747 0.3477984 1.8867719 3.2673090
#> [5,] 1.2791025 0.1161520 1.0705573 1.5282725
#> [6,] 0.7001766 0.1132798 0.5099103 0.9614382
max(abs(fitted(fit)[, "Estimate"] -
          frm_linpred(fit, type = "response")))
#> [1] 0