This a simple single-level Bernoulli model with weakly regularizing priors. Nothing fancy. Yet somehow its lppd is higher than the LL of the corresponding frequentist model, as seen below:

```
> baymod <- brm(y ~ x1 + x2 + x3 + x4, family = bernoulli, prior = prior(normal(0,2.5))+prior(normal(0,0.3), class = b, coef = x2), data = anon, iter = 5000, seed = 2023, cores = 4, backend = "cmdstanr")
# (x2 has smaller units than the rest, hence the tighter prior)
> baymod_loo <- loo(baymod)
> baymod_loo$elpd_loo + baymod_loo$p_loo
[1] -1453.191
```

Now for the frequentist version, which should have the highest possible log likelihood given that it is fit by maximizing the log likelihood:

```
> logLik(glm(y ~ x1 + x2 + x3 + x4, family = binomial, data = anon))
'log Lik.' -1453.282 (df=6)
```

How on earth is this possible? Beating a ML model at its own game, despite regularization by priors? The mind boggles. Can anyone help unboggle it? The anonymized datafile is attached.

anon.txt (29.0 KB)