Hi all,
I’ve fitted a multilevel model in brms and I want to know the MAP parameter values. My purpose is sampling from the MAP random effects structure to get the best guess of the predictive distribution without any estimation uncertainty. There is already great functionality to simulated new datasets from the random effects from randomly drawn posterior samples using brms::predict
and sample_new_levels = "gaussian"
. But, as I understand it, this incorporates the model uncertainty baked into the posterior samples. I’d like to do this for the MAP parameters.
I’ve already tried a couple of ways.

taking the MAP estimates of posterior samples using
bayestestR::map_estimate
. But when I simulate the predictive distribution by constructing the covariance matrix from the MAP sds and correlations, the resulting covariance matrix is not positive definite. I’m wonder whether this is because taking the MAP of each parameter value means that the combined parameters are no longer jointlycredible. 
extracting the stanmodel from brms using
mod[["fit"]]@stanmodel
and usingrstan::optimizing
. But this throws an error about variables not existing.
I’ve also wondered whether I’m overcomplicating things and it could be as simple as taking the posterior sample with maximum lp_
(if I interpret correctly that this would be the sample where the parameters match the model the most)?
Any direction appreciated.
Thanks for your time.
Callum