Having read that modeling measurement error in `brms`

via `me()`

terms are soft-deprecated in favor of `mi()`

, I’ve been trying to understand how each of these terms work in practice. In doing so, I’ve stumbled across some unexpected behavior (at least to me) of `emmeans`

(specifically `emtrends()`

) when using `mi()`

terms.

The issue arises when trying to calculate the effect of some predictor variable specified via `mi()`

using `emtrends()`

. Some observations are that (1) the result of `emtrends()`

changes every time, (2) the estimated trend deviates from the trend implied by the coefficient(s), and (3) the lower/upper HPD have an unexpectedly large range and magnitude. This does not seem to be an issue when the predictor is specified via `me()`

. Does someone have an explanation for this behavior?

Here’s a reproducible example:

```
library(brms)
library(emmeans)
set.seed(359)
# Simulate data
n <- 100
a <- 0
b <- 1
sigma <- 1
x <- rnorm(n, 0, 1)
xsd <- abs(rnorm(n, 0, 1))
xobs <- rnorm(n, x, xsd)
mu <- a + b*x
y <- rnorm(n, mu, 1)
df <- data.frame(y, xobs, xsd)
##### Model with me() #####
# Specify model formula
form_me <- bf(y ~ me(xobs, xsd))
# Fit model
fit_me <- brm(form_me, df, chains = 4, cores = 4)
# Model summary
summary(fit_me)
# emtrends (gives expected results)
# Specify arbitrarily small xsd to get effects without measurement error
emtrends(fit_me, ~ xobs, var = "xobs", at = list(xsd = 1e-10))
##### Model with mi() #####
# Specify model formula
form_mi <- bf(y ~ mi(xobs)) + bf(xobs | mi(xsd) ~ 1) + set_rescor(FALSE)
# Fit model
fit_mi <- brm(form_mi, df, chains = 4, cores = 4)
# Model summary
summary(fit_mi)
# emtrends (gives unexpected results)
emtrends(fit_mi, ~ xobs, var = "xobs", resp = "y", at = list(xsd = 1e-10))
```

`summary(fit_me)`

returns:

```
Family: gaussian
Links: mu = identity; sigma = identity
Formula: y ~ me(xobs, xsd)
Data: df (Number of observations: 100)
Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
total post-warmup draws = 4000
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept -0.03 0.12 -0.27 0.20 1.00 2620 2714
mexobsxsd 1.02 0.12 0.77 1.26 1.00 2407 2529
Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma 1.00 0.09 0.83 1.19 1.00 1965 2570
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).
```

`summary(fit_mi)`

returns:

```
Family: MV(gaussian, gaussian)
Links: mu = identity; sigma = identity
mu = identity; sigma = identity
Formula: y ~ mi(xobs)
xobs | mi(xsd) ~ 1
Data: df (Number of observations: 100)
Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
total post-warmup draws = 4000
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
y_Intercept -0.03 0.12 -0.26 0.20 1.00 5834 3015
xobs_Intercept 0.07 0.12 -0.17 0.31 1.00 5731 3294
y_mixobs 1.03 0.12 0.79 1.27 1.00 2903 3079
Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma_y 1.00 0.09 0.84 1.20 1.00 3787 2917
sigma_xobs 1.05 0.09 0.89 1.26 1.00 3657 3437
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).
```

The estimates from the two models are largely equivalent (except for the differences inherent to the different model specifications). However, using `emtrends()`

with the two models seems to produce different results:

`emtrends(fit_me, ~ xobs, var = "xobs", at = list(xsd = 1e-10))`

returns expected estimates:

```
xobs xobs.trend lower.HPD upper.HPD
0.0223 1.02 0.774 1.26
Point estimate displayed: median
HPD interval probability: 0.95
```

`emtrends(fit_mi, ~ xobs, var = "xobs", resp = "y", at = list(xsd = 1e-10))`

returns unexpected estimates:

```
xobs xobs.trend lower.HPD upper.HPD
0.0223 15 -1107 1202
Point estimate displayed: median
HPD interval probability: 0.95
```

System info:

- Operating System: macOS Ventura 13.4
- brms Version: 2.20.4
- emmeans Version: 1.10.0