# in- and out-of-sample predictions with mlogit using random parameters

I am using the mlogit R package to fit a mixed multinomial logit model -- that is, a multinomial logit model with random coefficients.

I fit my model with in-sample choice data consisting of in-sample individuals/decision makers. My model has a single random parameter which I have specified to be normally distributed. I obtain estimates for the mean and standard deviation of the random parameter's distribution. Now, if I call the predict() function on my mlogit object, how does the package compute predicted probabilities? Does it used individual-level estimates of the coefficients, or does it provide a Monte Carlo estimate with the estimated mean/variance of the random coefficient distribution?

Then, using the predict() function, I can also make predictions on new out-of-sample data with a new set of individuals. However, I can't find any documentation about how the package makes these predictions either.

I tried looking at the code, but it difficult to follow.

Here's a guess: using the estimated mean and standard deviation of the random parameter's distribution, the package might draw R coefficient samples from this fixed distribution and average the corresponding R logit probabilities. Does anyone know for sure if this is correct? Does it do this for both in- and out-of-sample estimates? If you are familiar with the code, please point me to the relevant file.

For example, if we assume one of the coefficient for decision-maker $h$ is normally distributed so that $\beta_h \sim \text{N}(\mu, \sigma^2)$, then calling mlogit() will produce estimates $\hat{\mu}$ and $\hat{\sigma}^2$.
With those estimates fixed, predictions are obtained through sampling to estimate the following integral representing the probability that individual $h$ chooses item $j$: $$p_{hj} = \int L(\beta_h) f(\beta_h | \hat{\mu}, \hat{\sigma}) d\beta_h$$ where $L(\beta_h)$ is the logit probability expressed as a function of $\beta_h$.