Hi,

I am trying to compute a distribution by integrating over all non-increasing categorical distributions of a given size $n$.

For instance, for $n=2$, each categorical distribution must follow $u_0 \geq u_1 \geq 0$ and $u_0 + u_1 = 1$. For $n=3$, it must follow $u_0 \geq u_1 \geq u_2 \geq 0$ and $u_0 + u_1 + u_2 = 1$.

For a given $n$, I want to integrate over all of the categorical distributions that are consistent with these constraints.

For example, if $n=2$, then I want to compute $p = \int_{U:~u_0 \geq u_1, u_0 + u_1 = 1} U \partial U$.

I was able to simplify this a bit by transforming:

$u_0 \in [\frac{1}{n}, 1]$

$u_1 \in [\frac{1-u_0}{n}, u_1]$

$u_2 \in [\frac{1-u_0-u_1}{n}, u_2]$

$\vdots$

$u_{n-1} = 1 - u_0 - u_1 - \cdots - u_{n-2}$

Using this, I can compute what I want using nested integrals:

$p_0 = \int_{1/2}^1 u_0 \partial u_0 = 3/8$

$p_1 = \int_{1/2}^1 1-u_0 \partial u_0 = 1/8$

If $n=3$, you'd compute

$p_0 = \int_{1/3}^1 \int_{(1-u_0)/2}^{u_0} u_0 \partial u_1 \partial u_0$

$p_1 = \int_{1/3}^1 \int_{(1-u_0)/2}^{u_0} u_1 \partial u_1 \partial u_0$

$p_2 = \int_{1/3}^1 \int_{(1-u_0)/2}^{u_0} 1-u_0-u_1 \partial u_1 \partial u_0$

I want to compute a handful of the $p_k$ $p_{k_1}, p{k_2}, \cdots $ for large $n$ (on the order of a few thousand). To be exact, all I *really* need are values proportional to these $p_k$.

I've managed to make an algorithm that can do this integration and compute a single $p_k$ in $O(n^3)$, but the constant isn't very fast, and it gets quite slow for large $n$.

This really feels like it *must* be a classic problem of some sort, but I just can't seem to see or find a way to get the closed form solution or a faster algorithm.

Also, I only need to do this numerically, **so sampling is an option**. Unfortunately, I could only work out a rejection sampling strategy (you throw away the entire configuration U as soon as one variable violates an inequality), so the number of iterations to generate one sample becomes roughly exponential in $n$. If you know of a way that I can sample nondecreasing uniform values without rejection, I'd love to know (I also tried generating $n$ uniform values and sorting them, but I don't know how to force them to sum to 1 and still keep the same distribution; scaling didn't seem to give me the right answer).

If this looks familiar to anything you've seen, I would really, really appreciate the help. I've been stuck on this for a few days, and I'm sure you can appreciate how that makes your brain itch like there's a little splinter you can't quite get out.

Oliver