I have completed a hierarchical multiple imputation of my data using the hmi package and have transformed the 10 imputed datasets it into a mitml.list. I am able to fit frequentists regression models to it using the with() function.

e.g., fittedn.model <- with(data, lmer(score ~ 1 + gender + (1|course_id:student)))

I can separate out one of the 10 imputed datasets and perform a bayesian regression. I could perform 10 independent regressions and then combine them and pool their variances using Rubins Rules, but it seems like there has to be a package that allows for the analysis of MI datasets. Can anyone point me in the right direction? Thanks!