# MCMC sampling for a model with a multinomial choice--so the parameters need to sum to 1

this is a head-scratcher for me, but a very interesting problem. So I have a stochastic simulation model for a hiring process. Basically different groups get hired into a company with different probabilities. So there is a multinomial process for hiring, where $$n$$ number of available positions go to members from $$m$$ number of groups in each year. So in a Bayesian sense the model would look like below--I left out the details for the prior parameters since they are not really important:

$$n \sim Binomial(\cdot, \cdot) \\ rate1 \sim TruncatedNormal(\cdot, \cdot)\\ rate2 \sim TruncatedNormal(\cdot, \cdot)\\ rate3 \sim TruncatedNormal(\cdot, \cdot)\\ rate4 \sim TruncatedNormal(\cdot, \cdot)\\ \\ \\ hires \sim Multinomial(n, [rate1, rate2, rate3, rate4]) \\$$

I have hiring data, and I was going to use Bayesian MCMC sampling to estimate the hiring rates for each group from the time-series data.

The challenge is that when I draw samples from the prior for each of the hiring rates, those rates have to sum to 1--otherwise the multinomial draw will not work. So I need to explore the set of hiring rates in such a way that all of the rates sum to 1.

I am not sure how to handle this kind of constraint in an MCMC scheme. Is there a distribution that I would sample the rates from such that they all sum to 1? It seems like a very interesting problem, but I only discovered it while debugging my own code, so I have not had much chance to think about it yet.

Any thoughts would be appreciated.

• Have you considered looking at the Dirichlet distribution? It is a probability distribution over discrete probability vectors. This would have the added bonus of Dirichlet-multinomial conjugacy. It might be important to check if that is the multinomial you had in mind, as it can be unclear whether authors are referring to the categorical or multinomial distribution. Jun 3, 2021 at 23:32
• I have extremely limited experience with MCMC, however, I have seen a paper before where the authors implement a sum-to-zero constraint trick in their scheme. Even though it's not sum-to-one, I will post it here when I have managed to find it, as it also shows how to implement this in code, and may be of interest. Jun 3, 2021 at 23:38
• @microhaus thanks for the suggestion here. I can think about the Dirichlet distribution, as that is a distribution over multinomials. Someone must have implemented something like that before, right. The other thing I can kind of think of is those convolutions in probability distributions, where you have the Probability of P(Z) given X + Z = 1. I think they solved those with characteristic functions or such, but it has been a while since I looked at that stuff :). Jun 4, 2021 at 0:08
• @microhaus you know, I did a little checking on the applications of the Dirichlet and it does seem to be used commonly for these types of problems. I have run into the Dirichlet before, but never had to use it myself. So that was a good call. Do you think it makes sense to keep this question open, or should I cancel it? Jun 4, 2021 at 1:22
• What I was referring to in my earlier comment is a blog post here. The blog uses the strategy from Baio, G., Blangiardo, M. (2011). Bayesian hierarchical model for the prediction of football results on a different data set. There are examples of how to code a sum-to-zero constraint (but not sum-to-one constraint for probability vectors) in both. IMO there is no need to close the question, you might receive concrete advice from MCMC experts. Jun 4, 2021 at 1:31

The problem does not seem to stand with MCMC but with the prior modelling. If the data comes from a Multinomial distribution $$\mathcal D_k(n,p_1,\ldots,p_k)$$ where the probability vector $$\mathbf{p}=(p_1,\ldots,p_k)$$ belongs to the $$k$$-dimensional simplex, the prior on $$\mathbf{p}$$ must put some mass on the simplex (and should logically be allocating all its mass to the simplex. Using truncated Normals is thus inadequate as the prior puts zero mass on the simplex.
As suggested by microhaus, a natural family of distributions in this setting is the family of Dirichlet priors. An alternative (or a generalisation) is to over-parametrise the model and set the constraint later. For instance, take $$p_1=\dfrac{\rho_1}{\sum_{i=1}^k \rho_i},\ldots,p_k=\dfrac{\rho_k}{\sum_{i=1}^k \rho_i}$$ with$$\rho_i\sim\mathcal N^+(a_i,b_i^2)$$ (when using a Gamma prior with a constant scale instead of a truncated Normal, this returns a Dirichlet).
• Since the prior is proper, the posterior is proper as well, hence MCMC has no convergence issue in theory. In the opposite, since the scale of the $\rho_i$'s is unconstrained by the data, the chain should converge better. Jun 4, 2021 at 9:03
• Yep, so I would have to use the normalized $p_k$s, as they are the only ones that fall within the probability simplex, I imagine. In other words, the $\rho_k$ seem to be thrown away, or just a means to get the normalized $p_k$ in the simplex. And that makes sense. Jun 4, 2021 at 14:58