Probability density function for the polar sine of uniformly distributed points on the sphere

Question

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?

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, but my attempts to generalise to $k=3$ do not match empirical distributions for $n=10$.

Answer 1

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.