Let
- $a = m -1$, $b = n - m - 1$, $d = a + b = n - 2$.
- $\theta_{ij}$, $1 \le i \ne j \le n$ be the angular separation between point $x_i$ and $x_j$.
- $\ell_m$ be the expected angular separation among a pair of $m^{th}$ nearest neighbors.
- $\mathcal{E}_m$ be the event that $x_2$ is the $m^{th}$ nearest neighbor of
$x_1$.
- $\mathbf{1}_m$ be the indicator function for event $\mathcal{E}_m$.
Since all points are equal, we can use any pair of points, say $x_1$ and $x_2$, as probe and $\ell_m$ will be equal to the expected value of $\theta_{12}$
subject to the constraint that event $\mathcal{E}_m$ happens. i.e.
$$\ell_m = \mathbf{E}( \theta_{12} | \mathcal{E}_m )
= \frac{\mathbf{E} ( \theta_{12}\mathbf{1}_m )}{\mathbf{P}( \mathcal{E}_m )}
= \frac{\mathbf{E} ( \theta_{12}\mathbf{1}_m )}{\mathbf{E}( \mathbf{1}_m )}
$$
For any $\theta \in [0,\pi]$, let $$p = \frac{1-\cos\theta}{2} \iff \theta = 2\sin^{-1}(p^{1/2})$$
When $\theta_{12} = \theta$, we have
$$\mathbf{P}( \theta_{1k} < \theta | \theta_{12} = \theta ) = \frac{1-\cos\theta}{2} = p
\quad\text{ for any } 2 \le k \le n
$$
Since these $n-2$ angular separations $\theta_{1k}$ are independent and there are
$\displaystyle\;\binom{d}{a}$ ways of picking $a$ out of $d$ points. we have
$$\mathbf{E}( \mathbf{1}_m | \theta_{12} = \theta ) = \binom{d}{a}p^a(1-p)^b$$
Since $x_2$ are distributed uniformly over the sphere, the probability for $\theta \le \theta_{12} \le \theta + d\theta$ is proportional to $\sin\theta_{12} d\theta_{12} \propto dp$. We have
$$
\mathbf{E}( \theta_{12}\mathbf{1}_m ) = \displaystyle\;\binom{d}{a}\int_0^1 p^a (1-p)^b \theta_{12} dp
\quad\text{ and }\quad
\mathbf{E}( \mathbf{1}_m ) = \displaystyle\;\binom{d}{a}\int_0^1 p^a (1-p)^b dp
$$
As a result,
$$\begin{align}
\ell_m
&= m\binom{n-1}{m}\int_0^1 p^a (1-p)^b 2\sin^{-1}(p^{1/2})\,dp\\
&= m\binom{n-1}{m}\sum_{k=0}^a(-1)^k \binom{a}{k} \int_0^1 (1-p)^{b+k} 2\sin^{-1}(p^{1/2})\, dp\\
&= m\binom{n-1}{m}\sum_{k=0}^a(-1)^{a-k} \binom{a}{k} \int_0^1 (1-p)^{d-k} 2\sin^{-1}(p^{1/2})\, dp
\end{align}
$$
For any $c \in \mathbb{N}$, we have
$$\begin{align}
\int_0^1 (1-p)^c 2\sin^{-1}(p^{1/2}) dp
&= \frac{-1}{c+1}\int_0^1 2\sin^{-1}(p^{1/2})\, d(1-p)^{c+1}\\
&= \frac{-1}{c+1}\left\{\bigg[ 2\sin^{-1}(p^{1/2}) (1-p)^{c+1} \bigg]_0^1 - \int_0^1 \frac{(1-p)^{c+1}}{\sqrt{p(1-p)}} dp \right\}\\
&= \frac{1}{c+1}\frac{\Gamma(\frac12)\Gamma(c+\frac32)}{\Gamma(c+2)}
= \frac{\pi}{(c+1)2^{2c+1}}\binom{2c+1}{c}
\end{align}
$$
Form this, we get
$$
\ell_m
= \pi m\binom{n-1}{m}\sum_{k=0}^{m-1}
\frac{(-1)^{m-1-k}}{(n-k-1)2^{2n-2k-3}}\binom{m-1}{k}
\binom{2n-2k-3}{n-k-2}
$$
In particular, the expected nearest neighbor angular separation $\displaystyle\;\ell_1 = \frac{\pi}{2^{2n-3}}\binom{2n-3}{n-2}$.