Projecting a distribution onto the unit sphere in an arbitrary norm

Question

So we represent the Dirichlet distribution as the projection of the $d$ independent gammas (on $R_+^d$ onto the unit simplex, and we arrive at that through the $L_1$ norm. That is, divide ${\bf x} \in R_+^{d}$ by $r = \sum x_i$ to push data onto the unit simplex, $S_{1}^{d-1}$. Then we integrate out $r$, and are left with Dirichlet (or generalized dirichlet, if using variable rate gammas).

I'm trying to set up the Jacobian of the transformation between ${\bf x} \in R_+^d$ and ${\bf y} \in S_{\omega}^{d-1}$, where $\omega$ specifies some arbitrary norm. If we want to narrow it down, then in the family of $p$-norms. That is, $$ r = \left(\sum x_i^p\right)^{1/p} $$ So you can imagine, $p = 1$ defines the simplex, $p = 2$ defines the hypersphere. My problem is specifically defined on the infinity norm; $\lim_{p\to\infty} (\sum x_i^p)^{1/p} = \max x_i$

The determinant of the Jacobian should still be $r^{d-1}$, but I'm not certain of that. It's really the ordering of variables that's messing me up; I need to end up with a random vector $(r,y_{(1)},\ldots,y_{(d-1)})$, and which dimension falls out is the max, but that's random.

Is it sufficient to show the determinant of the jacobian is the same no matter which dimension falls out?

EDIT: To be clear, the standard transformation of the is something like: $$ \begin{aligned} T(x_1,\ldots,x_d) &= (r, y_1,\ldots,y_{d-1})\\ \implies T^{-1}(r,{\bf y}) &= \left(ry_1,\ldots,ry_{d-1},r(1-\sum_{i=1}^{d-1}y_i^p)^{1/p}\right) \end{aligned} $$ But when I take the determinant of the Jacobian of that transformation, I end up with something like $$ r^{d-1}\left[(1 - \sum_{i = 1}^{d-1}y_i^p)^{1/p} + \sum_{i = 1}^{d-1}y_i^p / (1 - \sum_{l = 1}^{d-1}y_l^p)^{1/p - 1}\right]. $$ If I take $p\to\infty$, if anything but $y_d = 1$, then the denominator goes to 0.

Answer 1

Suppose $(X_1,\ldots,X_d)$ is a positive continuous random variable. Write its joint probability element as

$$f(x_1,\ldots,x_d)\,\mathrm{d}x_1 \cdots \mathrm{d}x_{d-1} \mathrm{d}x_d$$

so that the density function is $f(x_1,\ldots, x_d).$

Leaving the case of the $L^\infty$ norm for later (it is special), consider $0\le p\lt \infty.$ You ask about the density of the random variable

$$(Y_1,\ldots,Y_{d-1},R) = (X_1/R, \ldots, X_{d-1}/R, (X_1^p+\cdots+X_d^p)^{1/p}).$$

This is an invertible transformation because $X_i=R\,Y_i$ is defined for $i=1,2,\ldots,d-1$ and $X_d = R\left(1 - Y_1^p-\cdots-Y_{d-1}^p\right)^{1/p}$ also is unique.

Emulating the derivation of the Dirichlet distribution using the algebra of differential forms, compute

$$\mathrm{d}x_i = r\,\mathrm{d}y_i + y_i\,\mathrm{d}r,\ i=1,2,\ldots, d-1$$

and

$$\begin{aligned} px_d^{p-1}\,\mathrm{d}x_d &= d\left(x_d^p\right) = \mathrm d\left(r^p-x_1^p-\cdots-x_{d-1}^p\right)\\ &= p\left(r^{p-1}\,\mathrm{d}r - x_1^{p-1}\,\mathrm{d}x_1 - \cdots - x_{d-1}^{p-1}\,\mathrm{d}x_{d-1}\right). \end{aligned}$$

Plug these in to compute the old differential in terms of the new variables and their differentials:

$$\begin{aligned} \mathrm{d}x_1\wedge \cdots \wedge \mathrm{d}x_d &= \left(px_d^{p-1}\right)^{-1}\mathrm{d}x_1\wedge \cdots \wedge \mathrm{d}x_{d-1}\wedge \mathrm{d}(x_d^p)\\ &= \left(x_d^{p-1}\right)^{-1}\mathrm{d}x_1\wedge \cdots \wedge \mathrm{d}x_{d-1}\wedge \left(r^{p-1}\,\mathrm{d}r - x_1^{p-1}\,\mathrm{d}x_1 - \cdots - x_{d-1}^{p-1}\,\mathrm{d}x_{d-1}\right)\\ &= \left(x_d^{p-1}\right)^{-1}\mathrm{d}x_1\wedge \cdots \wedge \mathrm{d}x_{d-1}\wedge r^{p-1}\,\mathrm{d}r \\ &= \left(x_d^{p-1}\right)^{-1} \left(r\,\mathrm{d}y_1 + y_1\,\mathrm{d}r\right) \wedge \cdots \wedge \left(r\,\mathrm{d}y_{d-1} + y_{d-1}\,\mathrm{d}r\right) \wedge r^{p-1}\,\mathrm{d}r\\ &= r^{d-1}\left(\frac{r^{p-1}}{x_d^{p-1}}\right)\,\mathrm{d}y_1\wedge \cdots \wedge \mathrm{d}y_{d-1} \wedge \mathrm{d}r\\ &= r^{d-1}\left(1 - y_1^p - \cdots - y_{d-1}^p\right)^{1/p-1}\,\mathrm{d}y_1\wedge \cdots \wedge \mathrm{d}y_{d-1} \wedge \mathrm{d}r. \end{aligned} $$

The coefficient of the differential form is the Jacobian and (therefore) the density of $(Y_1,\ldots, Y_{d-1})$ is

$$\begin{aligned} &g(y_1,\ldots, y_{d-1}) \\ &= \int_0^\infty f\left(ry_1, \ldots, ry_{d-1}, r\left(1-y_1^p-\cdots-y_{d-1}^p\right)^{1/p}\right)\, r^{d-1}\left(1 - y_1^p - \cdots - y_{d-1}^p\right)^{1/p-1}\,\mathrm{d}r. \end{aligned}$$

---

The difficulty when $p=\infty,$ where $r=\max(x_i),$ is that $(y_1,\ldots, y_{d-1}) = (x_1/r,\ldots, x_{d-1}/r)$ are no longer coordinates for the image of the projection. These coordinates work only in the subset where $x_d$ is the largest.

A simple approach is to decompose the distribution according to the events determined by which of the $X_i$ is the largest. Because the joint distribution is continuous, the chance of a tie for largest is zero, so may neglect that possibility. Let $K$ be the random variable given by the index of the largest of the $X_i.$

Conditioning first on $K=d,$ we have $R=\max(X_i)=X_d$ and $Y_i = X_i/R = X_i/X_d$ for $i=1,2,\ldots, d-1.$ Computing as before we (easily) obtain

$$\mathrm{d}x_1\wedge \cdots \wedge \mathrm{d}x_{d} = r^{d-1}\, \mathrm{d}y_1\wedge \cdots \wedge \mathrm{d}y_{d-1}\wedge \mathrm{d}r.$$

Comparing this to the $L^p$ solution obtained previously, we may offer an intuitive (hand-waving) argument for this formula (as suggested in the question): as $p$ grows large, the factor $(1-y_1^p-\cdots -y_{d-1}^p)^{1/p-1}$ approaches $x_d^{-1} = r^{-1},$ giving $r^d(r^{-1})= r^{d-1}$ for the Jacobian.

This analysis applies with almost no change when conditioning on other values of $K:$ you only have to take the proper sign when integrating over $r$ in each case to assure the result is a positive density. In every case the Jacobian still equals $r^{d-1}.$ The full distribution of the projection is an equal mixture of the $d$ conditional distributions. Each of the $d$ components is supported on a different face of the unit cube $I_d = (0,1]^d\subset \mathbb{R}^d;$ the component corresponding to $K=k$ is supported on $\partial_k^{+} I_d = \{(y_1,\ldots,y_d)\mid 0\lt y_i\le 1;\ y_k=1\} \subset I_d.$