Why not sample a vector of size n, instead of sampling 4000 vectors of size 1 in a for loop?

This does not depend on i. You can do it outside of the for loop only once, instead of repeating the calculation 4000 times.

This is confusing. How are these probabilities? you mean parameters? but the parameters beta in the text were different: -2, 0, 2, 4.

Same as LP <- as.list(bp)? But why do you need the vector bp expressed as a list?

Not used subsequently.

I never used it, but I understand that you can use the categorical family of {brms} to fit multinomial models.

In the frequentist world, the package {mclogit} seems to fit this kind of models, and also the function multinom() from package {nnet}.

Hi! nice post, thank you! This is an interesting topic on which I worked a bit in the past trying to implement multinomial models in INLA using this trick. But I got lost in the details and computational problems and never got around to complete it.

And here is how to highlight and annotate text in the page. You can also add general comments for the full page.