Pooled post-processing of a frm_multiple() fit

brms's brm_multiple() combines the posterior draws of the m imputations, so fixef(), summary(), conditional_effects() and the other methods of its result summarize the combined draws. A frm_multiple() result is m maximum-likelihood fits, and these methods pool them by Rubin's rules with the Barnard-Rubin degrees of freedom instead:

Usage

# S3 method for class 'frmtmb_multiple'
fixef(
  object,
  summary = TRUE,
  robust = FALSE,
  probs = c(0.025, 0.975),
  pars = NULL,
  ...
)

# S3 method for class 'frmtmb_multiple'
summary(
  object,
  priors = FALSE,
  prob = 0.95,
  robust = FALSE,
  mc_se = FALSE,
  ...,
  waic = FALSE
)

# S3 method for class 'frmtmb_multiple'
conditional_effects(x, effects = NULL, ..., prob = 0.95)

# S3 method for class 'frmtmb_multiple'
vcov(object, correlation = FALSE, pars = NULL, ...)

# S3 method for class 'frmtmb_multiple'
coef(object, summary = TRUE, robust = FALSE, probs = c(0.025, 0.975), ...)

# S3 method for class 'frmtmb_multiple'
ranef(
  object,
  summary = TRUE,
  robust = FALSE,
  probs = c(0.025, 0.975),
  pars = NULL,
  groups = NULL,
  ...,
  condVar = FALSE
)

# S3 method for class 'frmtmb_multiple'
VarCorr(
  x,
  sigma = 1,
  summary = TRUE,
  robust = FALSE,
  probs = c(0.025, 0.975),
  ...
)

# S3 method for class 'frmtmb_multiple'
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
)

# S3 method for class 'frmtmb_multiple'
nobs(object, ...)

# S3 method for class 'frmtmb_multiple'
formula(x, ...)

# S3 method for class 'frmtmb_multiple'
family(object, resp = NULL, ...)

Arguments

object, x: A frmtmb_multiple from frm_multiple().
summary, robust: brms's arguments, refused when not the default: there are no draws.
probs: Probabilities of the quantile columns of fixef(), fitted() and VarCorr().
pars: Row names to keep, in the order given, as in brms. For ranef(), the coefficients to keep, as in ranef().
priors: If TRUE, summary() also prints the prior of the first imputation's fit.
prob: The coverage of the summary's and the curves' intervals.
mc_se: brms's argument, refused when TRUE.
waic: Old brms printed WAIC in summary() when it was TRUE. Maximum-likelihood fits have no WAIC, so TRUE is refused by name (brms 2.23.0 ignores the argument).
effects, ...: For conditional_effects(), passed to conditional_effects() on each imputation's fit. Otherwise refused by name.
correlation: If TRUE, vcov() returns the pooled correlation matrix.
groups: For ranef(), the grouping factors to keep, as in ranef().
condVar: For ranef(): if TRUE, each matrix carries the pooled conditional standard deviations as its "condSD" attribute.
sigma: Ignored, as in brms; nlme's generic carries it.
newdata, re_formula, scale, resp, dpar, nlpar: For fitted(): as in fitted.frmtmb_fit(), given to every imputation's fit. resp also names the response whose family family() returns.
ndraws, draw_ids, sort: For fitted(): brms's draws arguments, as on one fit.
allow_new_levels, sample_new_levels: For fitted(): as in fitted.frmtmb_fit(), given to every imputation's fit.

Value

fixef(): brms's coefficient matrix. summary(): a summary.frmtmb_multiple object, printed in brms's layout. conditional_effects(): a conditional-effects object, which plot() draws. vcov(), coef(), ranef(), VarCorr() and fitted(): the shape the same method gives on one fit. nobs(): a number.

Details

fixef() pools fixef() of each imputation, so it has brms's rows, names and order, ordinal thresholds included. Estimate is the mean over the imputations, Est.Error the square root of the total variance (within plus 1 + 1/m times between), and the quantile columns a t interval on the pooled degrees of freedom.
summary() prints that table, the distributional parameters without a formula pooled on their link scale and transformed back, and the pooled variance components, with the degrees of freedom and the fraction of missing information (fmi) of each row where brms prints R-hat and effective sample sizes.
conditional_effects() computes the curves of every imputation on the grid of the first (brms's grid, which brm_multiple() reads from the first data set), pools each grid point on the scale its Wald band is symmetric on (the link, or the logit of an ordinal category probability), and transforms back. A category probability is the exception for its estimate: that is the mean of the imputations' probabilities, so the categories sum to 1, and the band is moved out to it where needed. Only the Wald band pools this way: method = "predict" and the bootstrap and profile bands are refused by name.
vcov() is the total covariance of the fixef() rows: the mean of each imputation's vcov() plus 1 + 1/m times the covariance of the estimates between imputations, which is mice::pool()'s.
fitted() pools fitted() of each imputation row by row. It pools on the link scale of the predictor, as conditional_effects() does: each estimate goes through the link and its standard error through the delta method. The estimate and a t interval on the pooled degrees of freedom come back through the inverse link, so the interval stays inside the response's range, and Est.Error is the pooled standard error carried back by the delta method. A category probability takes its interval from the logit scale, but its estimate and Est.Error from Rubin's rules on the probability scale, so the categories sum to 1; an end on the wrong side of the estimate is moved to it. These are the values of the pooled conditional_effects(). scale = "linear" pools the linear predictor itself. In sample, each imputation predicts its own rows, so an imputed covariate enters at the value of each imputation; brms answers from the first imputed data set alone and warns that it does. With newdata, every imputation predicts the same rows.
coef() and ranef() are the mean over the imputations of the coefficients and conditional modes of each level, which is Rubin's point estimate. ranef(condVar = TRUE) pools the conditional standard deviations as Rubin's total variance: the mean conditional variance plus 1 + 1/m times the variance of the modes.
VarCorr() pools each standard deviation on the log scale and each correlation on the Fisher z scale, the scales frm_multiple() pools its variance components on, and transforms back, with Est.Error by the delta method. A covariance pools on its own scale. A component at its boundary (a standard deviation of zero in one imputation) has no log, and its row is NA.
emmeans::emmeans() pools the estimated marginal means of each imputation as emmeans pools a mice::mira: at the covariate means over the imputations, with the Barnard-Rubin degrees of freedom on the asymptotic reference that emmeans uses for one fit.
nobs() is the number of observations, which every imputation shares, and formula() and family() are those of the model every imputation is fitted with.

predict(), residuals(), pp_check(), logLik() (and so AIC() and BIC()), loo(), waic(), update() and model.frame() have no Rubin's-rules meaning. Each is refused by name, with what to call instead: anova() pools a model comparison, and x$fits[[j]] is the fit of imputation j.

brms's draws arguments have no meaning on pooled fits and are refused by name, as on one fit.

Examples

set.seed(8)
n <- 80
x <- rnorm(n)
y <- rnorm(n, 1 + 0.5 * x, 1)
x[sample(n, 15)] <- NA
imps <- lapply(1:3, function(i) {
  xi <- x
  xi[is.na(xi)] <- sample(x[!is.na(x)], sum(is.na(xi)), TRUE)
  data.frame(y = y, x = xi)
})
fm <- frm_multiple(bf(y ~ x), family = gaussian(), data = imps)
fixef(fm)
#>            Estimate Est.Error        Q2.5     Q97.5
#> Intercept 0.9240484 0.1330324 0.658971231 1.1891256
#> x         0.3156069 0.1495728 0.004987834 0.6262259
summary(fm)
#>  Family: gaussian 
#>  Links: mu = identity; sigma = log
#> 
#> Formula: y ~ x 
#>    Data: the imputed data sets (Number of observations: 80) 
#>  Method: ML on each of 3 imputations, pooled by Rubin's rules (Barnard-Rubin df)
#> 
#> Regression Coefficients:
#>           Estimate Est.Error l-95% CI u-95% CI    df  fmi
#> Intercept     0.92      0.13     0.66     1.19 73.93 0.01
#> x             0.32      0.15     0.00     0.63 21.49 0.24
#> 
#> Further Distributional Parameters:
#>       Estimate Est.Error l-95% CI u-95% CI    df  fmi
#> sigma     1.17      0.09        1     1.38 70.37 0.03
# the total covariance, whose diagonal is the squared Est.Error
vcov(fm)
#>             Intercept           x
#> Intercept 0.017697614 0.003125799
#> x         0.003125799 0.022372008
fitted(fm, newdata = data.frame(x = c(-1, 0, 1)))
#>       Estimate Est.Error      Q2.5     Q97.5
#> [1,] 0.6084416 0.1838968 0.2387417 0.9781414
#> [2,] 0.9240484 0.1330324 0.6589712 1.1891256
#> [3,] 1.2396553 0.2152237 0.8024554 1.6768552
ce <- conditional_effects(fm, "x")
head(ce$x[, c("x", "estimate__", "lower__", "upper__")])
#>           x   estimate__    lower__   upper__
#> 1 -3.011052 -0.026260132 -0.9510917 0.8985715
#> 2 -2.960928 -0.010440855 -0.9203149 0.8994332
#> 3 -2.910805  0.005378423 -0.8895595 0.9003163
#> 4 -2.860682  0.021197700 -0.8588266 0.9012220
#> 5 -2.810558  0.037016977 -0.8281175 0.9021515
#> 6 -2.760435  0.052836255 -0.7974335 0.9031060

Usage

Arguments

Value

Details

See also

Examples