How to sample $n$ observations from a multinomial distribution using binomial (or poisson) sampling?

Question

Context

I have $n$ observations which I'd like to sample with replacement for the purpose of bootstrap. A way to think about it is that we have a multinomial distribution with $n$ classes and that we'd like to draw from it $n$ items.

However, due to computational context, it's easier for me to use either a binomial or poisson distribution. So an alternative (known) way to do approximate sampling will be to count how many times we wish to have each item, but drawing a Binomial($n$, $\frac{1}{n}$) frequency for each of the items. When $n$ is large, we can approximate this using Poisson(1) (i.e.: poisson bootstrap).

The question

I'd like to use these sampling to get exact multinomial distribution, instead of an approximate one. How can this be achieved using either binomial or poisson distribution?

The issue is that since we use, for each of the items, Binomial($n$, $\frac{1}{n}$), we have a non-zero chance that the sum of these binoms will be different than $n$. If it is smaller then n, then we can just use the values up to the $n-1$ item, and for the $n$th item just use $n$ minus the sum of the previous draws. But if the sum is larger than $n$, this method won't work.

An idea I had here is to either use some sort of rejection sampling, by which whenever the sample is not exactly $n$, to just try it again (which feels VERY expensive).

Another method I thought about was to randomly shuffle the order of the $n$ items, and then take the first $k$ items, so that their sum would be equal to $n$, and if the $n$ items would sum to a smaller number than $n$, then to pick the $n$th value and have it be n minus the sum of the values up to that point.

The issue is that I don't have a proof that this method would yield a valid multinomial distribution (I suspect it will, but no proof).

Any suggestions/references would be lovely. Thanks upfront!

Answer 1

You can do it by progressing conditionally through the categories. I'm going to work from the last category backward (for a particular reason) but it can be done in any order as long as you're consistent in how you go about it.

Equal probability case:

Sample the count in the $n$-th category, $X_n\sim\text{bin}(n,\frac{1}{n})$.

Sample the count in the $n-1$-th category conditionally on $X_n=x_n$, i.e. $X_{n-1}\sim\text{bin}(n-x_n,\frac{1}{n-1})$.

Sample the count in the $n-2$-th category conditionally on $X_{n-1}+X_n=s_{n-1:n}$, i.e. $X_{n-2}\sim\text{bin}(n-s_{n-1:n},\frac{1}{n-2})$, ... and so on. Naturally, if $s_{i:n}$ hits $n$, all lower $X$'s are $0$ for that sample.

Unequal probability case:

Sample the count in the $n$-th category, $X_n\sim\text{bin}(n,p_n)$.

Sample the count in the $n-1$-th category conditionally on $X_n=x_n$, i.e. compute $p_{n-1:n} = \frac{p_{n-1}}{{1-p_n}}$ (i.e. scaling up the remaining probabilities since we no longer have the last category) and draw $X_{n-1}\sim\text{bin}(n-x_n,p_{n-1:n})$.

Sample the count in the $n-2$-th category conditionally on $X_{n-1}+X_n=s_{n-1:n}$, i.e. compute $p_{n-2:n} = \frac{p_{n-2}}{{1-p_{n-1}-p_n}}$ and draw $X_{n-2}\sim\text{bin}(n-s_{n-1:n},p_{n-2:n})$, ... and so on. Naturally, if $s_{i:n}$ hits $n$, all lower $X$'s are $0$ for that sample.

So you progress through, adjusting both the "n" and "p" to account for the conditioning on what categories are already drawn and what's left to draw from.

If there are a very large number of categories and so typically very small values of $p_i$ you may need to pay attention to numerical error in the scaling of those $p$ values as you progress.