me(x, sdx) in a formula says that the column x is a noisy
measurement of a latent, noise-free value, with known measurement
standard deviation sdx, and that the model uses the latent value.
It is brms's me() term with brms's meaning.
Value
me() is a formula term, not a free-standing function:
bf() reads it at parse time, and the value it contributes is the
noise-free predictor's coefficient together with the latent values
and the measurement-model hyperparameters, reachable through
fixef() and the Noise-free Terms section of summary(). This
page itself documents the term grammar and returns nothing.
Details
The model is $$x_i \sim N(\tilde{x}_i, sdx_i), \qquad \tilde{x}_i \sim N(meanme, sdme)$$
where meanme and sdme are estimated. The coefficient of the term
is bsp_me<x><sdx> in brms's names, for example bsp_mexsx for
me(x, sx).
me() is not a function. It is a formula special that the parser
reads, like mo() and mi().
Arguments of me()
x: the noisy variable, a numeric column or an expression of columns, such aslog(x).sdx: the known measurement SD, a column, an expression or a positive constant.gr: optional. The name of a grouping column. There is then one latent value per level ofgr, andxandsdxmust be constant within each level.
Estimation
The latent values are integrated out by the Laplace approximation,
together with the random effects. For a gaussian response with a
linear me() term the marginal likelihood is multivariate normal,
and the approximation is exact. An interaction of two me() terms,
or a non-gaussian response, makes it approximate.
Several me() terms that share a grouping have correlated latent
values by default, as in brms. Use set_mecor() to make them
independent.
me() can appear on its own, in a : or * interaction with
numeric variables or with other me() terms, and in the formula of
any distributional or nonlinear parameter. It cannot appear in a
group-level term, in a nonlinear formula body, inside another
function call, or in the same interaction as mo() or mi().
Parameters and methods
The hyperparameters are the template components meanme, logsdme
(the log of sdme) and thetame (the correlation parameters). The
Noise-free Terms section of summary() reports them on brms's scale
as meanme_<coef>, sdme_<coef> and corme__<coef1>__<coef2>.
Class "b" addresses the me() coefficients by their names, for
example coef = "mexsx".
Priors
The hyperparameters have brms's prior classes. default_prior()
lists brms's rows for them:
"meanme": the latent mean, one row perme()variable under a class row (coef = "mex"forme(x, sdx))."sdme": the latent standard deviation itself, bounded below by 0, with the same rows. The density carries the log-Jacobian of the log map, as class"sd"does."corme":lkj(eta)on the correlation matrix of the terms that share agr =, one class row and, for a grouped set, a row withgroup =its grouping variable.
On frm() a prior on them makes the fit a MAP estimate, and the
defaults apply nothing. frmtmb.sample::frm_sample() applies brms's
defaults: flat on meanme and on sdme, and lkj(1) on corme.
The three classes take no resp, dpar or nlpar, because one set
of latent values serves the whole model, as in brms.
anova() compares only fits with the same me() calls, because each
call puts its noisy variable's measurements into the likelihood.
AIC() does not check this, so compare by AIC only across fits that
share their me() calls.
fitted(), predict(), residuals() and simulate() on the fitted data
use the estimated latent values (the conditional modes). A new row has no
latent value, and brms draws one from N(x, sdx) per draw. So do predict()
and simulate(), per replicate. fitted() reports the mean of the expected
response over the draw, with its variance in Est.Error and its quantiles in
the Q columns. frm_linpred() stays at the observed value and adds the
variance to se.fit. New data must hold sdx, as in brms. Two responses or
two dpars read one latent value, a stated divergence from brms, which draws
one per predictor. Thus a category probability moves, as it is a joint
function of both predictors: 0.3336 against brms's 0.3106. A row with no x
or sdx reads NA, with a warning of class frmtmb_me_missing. brms reads
NaN and warns.
See also
set_mecor(); mi() for a measurement model that has
covariates of its own, bf(x | mi(sdx) ~ z).
Examples
set.seed(1)
n <- 200
tx <- rnorm(n, 1, 0.8)
d <- data.frame(x = tx + rnorm(n, 0, 0.4), sx = 0.4)
d$y <- 2 + 0.7 * tx + rnorm(n, 0, 0.5)
fit <- frm(bf(y ~ me(x, sx)) + gaussian(), data = d)
fixef(fit)
#> Estimate Est.Error Q2.5 Q97.5
#> Intercept 1.9213953 0.08471423 1.7553585 2.0874322
#> mexsx 0.7443322 0.06947982 0.6081543 0.8805102
# the naive slope is attenuated by the measurement error
coef(lm(y ~ x, data = d))[["x"]]
#> [1] 0.5741128
# brms's prior rows for the latent mean and sd
default_prior(bf(y ~ me(x, sx)) + gaussian(), data = d)
#> route = "fit": the prior defaults frm() applies
#> prior class coef group resp dpar nlpar lb ub source
#> 1 (flat) Intercept NA NA default
#> 2 (flat) b NA NA default
#> 3 (flat) b mexsx NA NA default
#> 4 (flat) meanme NA NA default
#> 5 (flat) meanme mex NA NA default
#> 6 (flat) sdme 0 NA default
#> 7 (flat) sdme mex NA NA default
#> 8 (flat) sigma NA NA default
fitp <- frm(bf(y ~ me(x, sx)) + gaussian(), data = d,
prior = set_prior("normal(1, 1)", class = "meanme") +
set_prior("exponential(1)", class = "sdme"))
prior_summary(fitp)
#> prior class coef group resp dpar nlpar lb ub source
#> normal(1, 1) meanme <NA> <NA> user
#> exponential(1) sdme <NA> <NA> user