Let $X_1$ and $X_2$ be uniform on $n$-spheres. What is the distribution of $\| X_1+X_2\|$?

Solution 1:

Given $U=X_1+X_2$ in $\mathbb{R}^n$ where $X_i$ are random points on the $n-1$-spheres $||X_i||=r_i$, and $R=||U||$, we have $$ R^2 = U^2 = X_1^2 + X_2^2 + 2X_1\cdot X_2 = r_1^2 + r_2^2 + 2r_1r_2\cos\Theta $$ where $\Theta\in[0,\pi]$ is the angle between $X_1$ and $X_2$. So, $\Theta$ is a random variable corresponding to the angle between two random points on the $n-1$-sphere.

Let's start tackling $\Theta$ and $\cos\Theta$ directly. Note first that these do not depend on the lengths $||X_i||=r_i$. Also, if we pick $X_1$ first, we can either rotate or choose a coordinate system so $X_1=[1,0,\ldots,0]$: ie, $X_1$ points along the first axis (aka $x$-axis in low dimensions). This is true because $X_2$ is uniformly distributed and independent of $X_1$. So, basically, the distribution of $\Theta$ (or $\cos\Theta$) is the same as the angle between a random point on the unit $n-1$-sphere and $[1,0,\ldots,0]$.

Now, let $Z = [Z_1,\ldots,Z_n]$ be a random point on the unit $n-1$-sphere: ie, so that $Z_1^2+\cdots+Z_n^2=1$. Then, $\cos\Theta=Z\cdot[1,0,\ldots,0]=Z_1$. So what we are after is the distribution of $Z_1$ for random points $Z$ on the unit $n-1$-sphere.

We can express the $n-1$-dimensional area of the unit $n-1$-sphere as $$ \omega_{n-1} = \int_0^\pi \omega_{n-2}(\sin\theta)^{n-2}\,d\theta $$ where $\omega_{n-2}(\sin\theta)^{n-2}$ is the $n-2$-area of an $n-2$-sphere with radius $\sin\theta$. Since we are after a uniform probability distribution, we need to divide this by $\omega_{n-1}$.

Next, we wish to express this in terms of the coordinate $z_1=\cos\theta$, which, using $d\theta/dz_1=-\sin\theta$ and $\sin\theta=\sqrt{1-z_1^2}$, gives us $$ \begin{align} \int_0^\pi \frac{\omega_{n-2}}{\omega_{n-1}} (\sin\theta)^{n-2}\, d\theta &=\int_{-1}^1 \frac{\omega_{n-2}}{\omega_{n-1}} (\sin\theta)^{n-2}\, \left|\frac{d\theta}{dz_1}\right|\,dz_1 \\ &=\int_{-1}^1 \frac{\omega_{n-2}}{\omega_{n-1}} (1-z_1^2)^{\frac{n-3}{2}}\, dz_1. \end{align} $$ Replace the boundary $[-1,1]$ for $z_1$ with any other interval, and you get the probability of $Z_1=\cos\Theta$ within that interval; so the probability density of $Z_1=\cos\Theta$ is $$ f_{\cos\Theta}(z) = \frac{\omega_{n-2}}{\omega_{n-1}} (1-z^2)^{\frac{n-3}{2}}. $$ This is basically just stating that for random variables $Y=h(X)$, the probability densities are related by $f_X(x)=f_Y(y)\cdot\left|h'(x)\right|$.

Returning to $R$, we already know $R^2$ is linear in $\cos\Theta$ with values in $[(r_1-r_2)^2, (r_1+r_2)^2]$. Entering the distribution of $\cos\Theta$, this gives the density of $S=R^2$: $$ f_{R^2}(s) = \frac{\omega_{n-2}}{2r_1r_2\omega_{n-1}} \left[ 1 - \left(\frac{s-r_1^2-r_2^2}{2r_1r_2}\right)^2 \right]^{\frac{n-3}{2}}. $$ Now, $f_R(r) = 2rf_{R^2}(r^2)$ (same rule as above for change of variables) which yields $$ f_{R}(r) = \frac{r\omega_{n-2}}{r_1r_2\omega_{n-1}} \left[ 1 - \left(\frac{r^2-r_1^2-r_2^2}{2r_1r_2}\right)^2 \right]^{\frac{n-3}{2}}. $$

Note that for points uniformly distributed between two radii, you should have density $f_R(r)=ar^{n-1}$ for some constant $a$ as the $n-1$-areas of the $n-1$-spheres of radius $r$ is $\omega_{n-1}r^{n-1}$. So for no $n$ is that the case.

This gives the distribution $F_R(r)$ of the distance from the origin. The probability density of the $n$-dimensional vector $U$ is found by dividing by the $n-1$-area $\omega_{n-1}r^{n-1}$ of the $n-1$-sphere of radius $r=||u||$: $$ f_{U}(u) = \frac{F_R(||u||)}{\omega_{n-1} ||u||^{n-1}} = \frac{\omega_{n-2}}{||u||^{n-2}r_1r_2\omega^2_{n-1}} \left[ 1 - \left(\frac{||u||^2-r_1^2-r_2^2}{2r_1r_2}\right)^2 \right]^{\frac{n-3}{2}}. $$

As for the $k$-area of the unit $k$-sphere, $$ \omega_k = \frac{2\pi^{\frac{n+1}{2}}}{\Gamma\left(\frac{n+1}{2}\right)}, $$ as may be found on Wikipedia, where the gamma-function satisfies $\Gamma(s+1)=s\Gamma(s)$, and $\Gamma(n+1)=n!$ for integers.

As a side-note, the $1-z^2$ term inside the brackets may be rewritten $$ 1 - \left(\frac{r^2-r_1^2-r_2^2}{2r_1r_2}\right)^2 = \frac{[r^2-(r_1-r_2)^2]\cdot[(r_1+r_2)^2-r^2]}{(2r_1r_2)^2} $$ which helps highlight that $r^2$ lies between $(r_1-r_2)^2$ and $(r_1+r_2)^2$.

Solution 2:

It is known that the uniform distribution on the unit $(n-1)$-sphere can be represented as a standard multivariate Gaussian divided by its norm. Therefore, \begin{align*} X_1 \overset{\mathcal{D}}{=} r_1\frac{Z_1}{\|Z_1\|} \qquad \text{and} \qquad X_2 \overset{\mathcal{D}}{=} r_2\frac{Z_2}{\|Z_2\|} \end{align*} where $Z_1, Z_2 \overset{\text{iid}}{\sim} N(\mathbf{0}_n, I_{n\times n})$. Hence, \begin{align*} \|U\| \overset{\mathcal{D}}{=} \left\|r_1\frac{Z_1}{\|Z_1\|} + r_2\frac{Z_2}{\|Z_2\|}\right\| = \sqrt{r_1^2 + r_2^2 + 2r_1r_2 \frac{Z_1^\intercal Z_2}{\|Z_1\|\|Z_2\|}} \end{align*} The next step is to find the distribution of $P \overset{\text{def}}{=} \frac{Z_1^\intercal Z_2}{\|Z_1\|\|Z_2\|}$, which takes the form of Pearson's correlation coefficient, except we are not subtracting out the sample mean in the variance/covariance calculation. You can actually show that

\begin{align*} T \overset{\text{def}}{=} \frac{P}{\sqrt{1-P^2}} \sim t_{n-1}/\sqrt{n-1} \end{align*} where $t_n$ is the t-distribution with $n$ degrees of freedom. This follows from the proof in Hotelling's "New Light on the Correlation Coefficient and its Transforms" (1953), changed from $n-2$ to $n-1$ because of not needing to estimate the mean.

Solution 3:

Continuing on from Tom Chen's answer, let $X\sim t_{n-1}$ and $$ f(x) = \frac x{(1-x^2)^{1/2}}. $$ Then $f(P)=_d (n-1)^{-1/2}X$, so that $\|U\|=_d g(X)$, where $$ g(x) = \sqrt{r_1^2 + r_2^2 + 2r_1r_2 f^{-1}((n-1)^{-1/2}x)}\,. $$ Note that $X$ has density $$ \varphi(x) = C_n\left({1 + \frac{t^2}{n-1}}\right)^{-n/2}, $$ where $$ C_n = \frac{\Gamma(n/2)}{\sqrt{(n-1)\pi}\,\Gamma((n-1)/2)}. $$ Thus, the density of $\|U\|$ is $$ h(y) = \left({\frac d{dy}(g^{-1}(y))}\right)\varphi(g^{-1}(y)), $$ defined for $r_1-r_2\le y\le r_1+r_2$. If we let $$ \gamma_y = \frac{y^2 - r_1^2 - r_2^2}{2r_1r_2}, $$ then $g^{-1}(y)=\sqrt{n-1}\,f(\gamma_y)$, so that \begin{align} \frac d{dy}(g^{-1}(y)) &= \sqrt {n-1}\,f'(\gamma_y)\left({\frac y{r_1r_2}}\right)\\ &= \frac{\sqrt {n-1}}{r_1r_2}\frac y{(1 - \gamma_y^2)^{3/2}}. \end{align} Also, \begin{align} \varphi(g^{-1}(y)) &= C_n(1 + f^2(\gamma_y))^{-n/2}\\ &= C_n\left({\frac1{1 - \gamma_y^2}}\right)^{-n/2}\\ &= C_n(1 - \gamma_y^2)^{n/2}. \end{align} Thus, $$ h(y) = \frac{\sqrt {n-1}\,C_n}{r_1r_2}\,y\,(1 - \gamma_y^2)^{(n-3)/2}. $$ Putting it all together, the density of $\|U\|$ is $$ h(y) = \frac{\Gamma(n/2)}{r_1r_2\sqrt\pi\,\Gamma((n-1)/2)} \,y\,\left({1 - \left({ \frac{y^2 - r_1^2 - r_2^2}{2r_1r_2} }\right)^2}\right)^{(n-3)/2} $$ for $r_1-r_2\le y\le r_1+r_2$, and $0$ otherwise.

EDIT:

Regarding the comment, when $n=3$, since $\Gamma(3/2)=\sqrt\pi\,/2$ and $\Gamma(1)=1$, this reduces to $$ h(y) = \frac1{2r_1r_2}y. $$ We then have \begin{align} \int_{r_1-r_2}^{r_1+r_2} h(y)\,dy &= \frac1{4r_1r_2}y^2 \bigg|_{r_1-r_2}^{r_1+r_2}\\ &= \frac1{4r_1r_2}((r_1+r_2)^2 - (r_1-r_2)^2)\\ &= \frac1{4r_1r_2}(4r_1r_2) = 1. \end{align}