If I sample three points independently, uniformly at random on an $n$-dimensional sphere of radius $R$, what is the probability density function of their [polar sine](https://en.wikipedia.org/wiki/Polar_sine)?


More generally, for $k<n$ if I sample $k$ points independently uniformly at random on an $n$-dimensional sphere of radius $R$, what is the probability density function of their polar sine?


An asymptotic form for large $n$ can be deduced in the case $k=2$ from Lemma 4 of Appendix A of Laarhoven's [Sieving for shortest vectors in lattices using angular locality-sensitive hashing](https://eprint.iacr.org/2014/744.pdf), but my attempts to generalise to $k=3$ do not match empirical distributions for $n=10$.



Probability density function for the polar sine of uniformly distributed points on the sphere
For $k=3$, $0\le s\le 1$ it is
$$c\_ns^{n-3}\left(\frac\pi2-\arcsin(s)\right)$$
where
$$c\_n=\begin{cases}\frac2\pi(n-2)\frac{2^{n-2}\left(\frac{n-2}2\right)!}{(n-2)!}&\text{ for $n$ even}\\
\frac{(n-2)!}{2^{n-3}\left(\frac{n-3}2\right)!}&\text{ for $n$ odd}.\end{cases}$$


Using polyspherical coordinates, the cumulative density function for fixed a given sine $s\_0$ between the first two vectors is
$$\frac{A\_1A\_{n-3}}{\frac12B(1,(n-1)/2)A\_{n-1}}\int\_0^u\cos\theta\sin^{n-3}\theta$$
where
$$u=\begin{cases}\arcsin\left(\frac s{s\_0}\right)&\text{ for $s<s\_0$}\\
\frac\pi2&\text{ for $s\ge s\_0$}\end{cases}.$$


Integrating this over $s\_0=\sin\theta$ gives a cumulative distribution function for $s$ whose derivative is the pdf.


