# Marginal density of uniform spherical distribution

Suppose that $$X$$ is distributed uniformly in the scaled $$n$$-sphere $$\sqrt{n} \mathbf{S}^{n-1} \subset \mathbf{R}^n$$. Then apparently the distribution of $$(X_1, \dots, X_k)$$, the first $$k < n$$ coordinates of $$X$$ has density $$p(x_1, \dots, x_k)$$ with respect to Lebesgue measure in $$\mathbf{R}^k$$, moreover if $$r^2 = x_1^2 + \cdots + x_k^2$$, then it is proportional to $$\left(1 - \frac{r^2}{n}\right)^{(n-k)/2 - 1}, \quad \text{if}~0 \leq r^2 \leq n,$$ and otherwise is 0. I tried to compute this using the fact that $$(X_1, \dots, X_k) \stackrel{\rm d}{=} \sqrt{n} (g_1, \dots, g_k)/\sqrt{g_1^2 + \cdots + g_n^2}$$, when $$g_i$$ are iid standard normal variables, but it was somewhat unclear to me how to use even this representation to compute the density. Can anyone sketch the details for me?

$$\newcommand{\R}{\mathbb{R}} \newcommand{\x}{\mathbf{x}} \newcommand{\X}{\mathbf{X}}$$ This is to present a formalization of the answer by Carlo Beenakker, without explicit use of the delta function.

We are going to assume that $$\X=(X_1,\dots,X_n)$$ is uniformly distributed on the unit sphere $$\mathbb S^{n-1}$$, rather than on $$\sqrt n\,\mathbb S^{n-1}$$.

For each real $$t\in(0,1)$$, define the measures $$\mu_t$$ and $$\nu_t$$ over $$\R^n$$ by the conditions $$\begin{equation*} \int f\,d\mu_t=\int_{\R^n}d\x f(\x)1_{1-t<|\x|^2\le1} \end{equation*}$$ and $$\begin{equation*} \int f\,d\nu_t=\frac{\int f\,d\mu_t}{\int d\mu_t} \end{equation*}$$ for all (say) nonnegative continuous functions $$f\colon\R^n\to\R$$, where $$|\cdot|$$ denotes the Euclidean norm. Then $$\nu_t$$ is a probability measure converging (as $$t\downarrow0$$) to the Haar measure $$h$$ on the unit sphere in $$\R^n$$, in the sense that (say) $$\begin{equation*} \int f\,d\nu_t\to\int f\,dh \end{equation*}$$ for all nonnegative continuous functions $$f\colon\R^n\to\R$$.

Take now any function $$f\colon\R^n\to\R$$ such that $$\begin{equation*} f(\x)=g(\x_{n-1}) \end{equation*}$$ for some nonnegative continuous function $$g\colon\R^{n-1}\to\R$$ and all $$\x\in\R^n$$, where $$\x_j:=(x_1,\dots,x_j)$$ for $$\x=(x_1,\dots,x_n)\in\R^n$$ and $$j=1,\dots,n-1$$. Then
\begin{align*} \int f\,d\mu_t&=\int_{\R^{n-1}}d\x_{n-1}\,g(\x_{n-1})\int_\R du\, 1_{1-t<|\x_{n-1}|^2+u^2\le1} \\ &=\int_{\R^{n-1}}d\x_{n-1}\,g(\x_{n-1}) (1+o(1))t\,(1-|\x_{n-1}|^2)^{-1/2}\,1_{|\x_{n-1}|<1}. \end{align*} Also, $$\begin{equation*} \int d\mu_t=\int_{\R^n}d\x\, 1_{1-t<|\x|^2\le1}\propto(1+o(1))t, \end{equation*}$$ where $$\propto$$ means an equality up to a constant factor, depending only on $$n,k$$. So, $$\begin{equation*} \int f\,dh=\lim_{t\downarrow0} \int f\,d\nu_t\propto\int_{\R^{n-1}}d\x_{n-1}\,g(\x_{n-1}) (1-|\x_{n-1}|^2)^{-1/2}\,1_{|\x_{n-1}|<1}. \end{equation*}$$ Thus, the joint pdf of $$\X_{n-1}=(X_1,\dots,X_{n-1})$$ is given by $$\begin{equation*} p_{n-1}(\x_{n-1})\propto(1-|\x_{n-1}|^2)^{-1/2}\,1_{|\x_{n-1}|<1}. \end{equation*}$$ Now successively integrating $$p_{n-1}(\x_{n-1})$$ ($$n-1-k$$ times) in $$x_{n-1},\dots,x_{k+1}$$ and each time using the formula $$\begin{equation*} \int_0^{b^{1/2}}(b-u^2)^p du=c_p b^{p+1/2} \end{equation*}$$ for real $$b>0$$ and $$p>-1$$ with $$c_p:=\int_0^1(1-u^2)^pdu\in(0,\infty)$$ (so that $$1/2$$ is added to the exponent $$p$$ after such an integration), we see that the joint pdf of $$\X_k=(X_1,\dots,X_k)$$ is given by $$\begin{equation*} p_k(\x_k)\propto(1-|\x_k|^2)^{(n-k)/2-1}\,1_{|\x_k|<1}, \end{equation*}$$ as desired.

• This is a wonderfully simple calculation and was illuminating. Thank you for the follow up questions, your persistence, and both answers. Commented May 8, 2020 at 20:35
• @DrewBrady : Thank you for your kind words. However, I feel grateful to Carlo for teaching me something of value, with patience. Commented May 8, 2020 at 21:03
• Do you by chance know how the ratio$\frac{ \int^{b}_{a} (n-x^{2})^{(n-2)/2}dx}{\int^{\sqrt{n}}_{-\sqrt{n}} (n-x^{2})^{(n-2)/2}dx} = \frac{ \int^{b}_{a} (1-x^{2}/n)^{(n-2)/2}dx}{\int^{\sqrt{n}}_{-\sqrt{n}} (1-x^{2}/n)^{(n-2)/2}dx}$ is related to marginal density of the first coordinate? Commented Feb 13, 2022 at 15:55
• McKean in his paper: "Geometry of Differential Space" makes use of it to prove that the marginal density of the first coordinate tends to a standard Gaussian. The paper can be found here: projecteuclid.org/journals/annals-of-probability/volume-1/… Commented Feb 13, 2022 at 15:57
• @user135520 : Thank you for the reference to McKean. Yes. it is well known and easy to prove that the distribution of the square of a coordinate, say $U_1$, of a uniformly distributed unit random vector $U$ in $\mathbb R^n$ has the beta distribution with parameters $1/2,(n-1)/2$. The Gaussian approximation to the distribution of $U_1$ easily follows from the representation of $U$ in terms of $n$ iid standard Gaussians. Commented Feb 13, 2022 at 16:30

With $$X$$ uniformly distributed over the unit $$n$$-sphere, the joint probability distribution of all $$n$$ elements of $$X$$ is a Dirac delta function, $$P(X_1,X_2,\ldots X_n)\propto\delta\left(1-\sum_{j=1}^n X_j^2\right).\qquad\qquad(1)$$ Now you integrate out elements one by one, to obtain the marginal distribution $$P_k$$ of $$k elements. The first integration gives $$P_{n-1}(X_2,X_3,\ldots X_n)\propto\left(1-\sum_{j=2}^n X_j^2\right)^{-1/2}\theta\left(1-\sum_{j=2}^n X_j^2\right),\qquad(2)$$ with $$\theta$$ the unit step function. The second integration gives $$P_{n-2}(X_3,\ldots X_n)\propto\theta\left(1-\sum_{j=3}^n X_j^2\right),$$ the third integration $$P_{n-3}(X_4,\ldots X_n)\propto\left(1-\sum_{j=4}^n X_j^2\right)^{1/2}\theta\left(1-\sum_{j=4}^n X_j^2\right),$$ and so on. Each additional integration increases the power by 1/2, $$P_{k}(X_{n-k+1},\ldots X_n)\propto\left(1-\sum_{j=n-k+1}^n X_j^2\right)^{(n-k)/2-1}\theta\left(1-\sum_{j=n-k+1}^n X_j^2\right).$$ This is the answer in the OP (without rescaling the radius of the $$n$$-sphere, so $$r^2/n\mapsto r^2$$).

As requested in the comments, a more detailed exposition of the various steps. • **First step:** the delta function. Denote the surface measure on the unit $$n$$-sphere as $$d\Omega$$, and $$\int d\Omega=A_n$$ the surface area. Uniformity of a distribution on the unit $$n$$-sphere means uniformity with measure $$d\Omega$$. I maintain that the joint probability distribution of the components of the vector $${\mathbf X}=(X_1,X_2,\ldots X_n)$$, uniformly distributed on the unit $$n$$-sphere, is given by Eq. (1) with normalization constant $$2/A_n$$. Let us check this by calculating the expectation value of an arbitrary function $$f$$ of $${\mathbf X}$$. For that purpose I transform to hyperspherical coordinates $$r,\phi_1,\phi_2,\ldots\phi_{n-1}$$, $$\mathbb{E}[f(\mathbf{X})]=\int dX_1\int dX_2\cdots\int dX_n \,f(X_1,X_2,\ldots X_n)P(X_1,X_2,\ldots X_n)$$ $$\qquad=\int_0^\infty r^{n-1} dr \int d\Omega\, f(r,\phi_1,\phi_2,\ldots\phi_{n-1})\frac{2}{A_n}\delta(1-r^2)$$ $$\qquad=\frac{1}{A_n}\int d\Omega\, f(r=1,\phi_1,\phi_2,\ldots\phi_{n-1}).$$ In the last step I used that $$\int_0^\infty r^{n-1}\delta(1-r^2)\,dr=1/2$$ for $$n\geq 2$$. • **Second step:** integration of the delta function, to arrive at Eq. (2). From now on I will ignore the normalization constants, these can easily be recovered at the end. Let me abbreviate $$\sum_{j=2}^n X_j^2=s_2$$. The marginal distribution $$P_1(X_2,X_3,\ldots X_n)$$ is obtained by definition upon integration of $$P(X_1,X_2,X_3,\ldots X_n)$$ over $$X_1$$. I carry out this integration in cartesian coordinates, changing variables to $$q=X_1^2$$, $$P_1(X_2,X_3,\ldots X_n)\propto \int_{-\infty}^\infty dX_1\delta(1-s_2-X_1^2),$$ $$\qquad=\int_0^\infty\delta(1-s_2-q)\frac{dq}{\sqrt q}=(1-s_2)^{-1/2}\theta(1-s_2).$$ • **Third and following steps:** The following steps, subsequent integrations of $$X_2,X_3,\ldots$$ are now immediate consequences of the integral $$\int_0^a(a^2-x^2)^p\,dx=c_p a^{1+2p}.$$
• Hm. What does it mean that this distribution "is a delta function"? Commented May 7, 2020 at 16:24
• @ Carlo Beenakker - are you sure you understand what a joint distribution is?
– R W
Commented May 7, 2020 at 18:51
• certainly, let's work this out for $n=2$; we then have in polar coordinates $x_1=r\cos\phi$, $x_2=r\sin\phi$ and we would expect the uniform distribution on the unit circle to be $P(\phi)=\text{constant}$. So let's check that the delta function distribution indeed gives this: take $P(x_1,x_2)\propto\delta(1-x_1^2-x_2^2)$ and integrate to obtain $P(\phi)=\int_0^\infty P(r\cos\phi,r\sin\phi) rdr\propto\int_0^\infty \delta(1-r^2)rdr=\text{constant}$. The multivariate Gaussian is a convenient way to generate the random coordinates, but for this calculation the delta function is easier. Commented May 7, 2020 at 21:02
• Carlo, I'd really like to understand your solution, especially because it promises an improvement over the Gaussian approach. Alas, I don't understand almost any part of your answer. Especially hard for me to understand your step-by-step integration. Can you write down in detail the integrals for at least the first two steps and how you take them? Commented May 7, 2020 at 23:29
• I think I got it -- thank you. Conversations with physicists are not always easy, but usually fruitful. :-) Commented May 8, 2020 at 16:18

Let $$G_1,\dots,G_n$$ be iid standard normal random variables. Then the random vector $$\begin{equation*} (Y_1,\dots,Y_n):=\Big(\frac{G_1}{\sqrt{\sum_1^n G_j^2}},\dots, \frac{G_n}{\sqrt{\sum_1^n G_j^2}}\Big) \end{equation*}$$ is uniformly distributed on the unit sphere $$\mathbb S^{n-1}$$. Let \begin{equation*} \begin{aligned} Z_i&:=Y_i=\frac{G_i}{\sqrt{\sum_1^n G_j^2}}&\text{ if }i\le k,\\ Z_i&:=G_i&\text{ if }i> k. \end{aligned} \end{equation*} We want to find the joint pdf of $$(Z_1,\dots,Z_k)$$, which is the same as the joint pdf of $$(X_1,\dots,X_k)/\sqrt n$$.

The vector $$(Z_1,\dots,Z_n)$$ is obtained from $$(G_1,\dots,G_n)$$ by the transformation given by \begin{equation*} \begin{aligned} z_i&:=\frac{g_i}{\sqrt{\sum_1^n g_j^2}}&\text{ if }i\le k,\\ z_i&:=g_i&\text{ if }i> k. \end{aligned} \end{equation*} The transformation inverse to this is given by \begin{equation*} \begin{aligned} g_i&:=\sqrt{s_2}\frac{z_i}{\sqrt{1-s_1}}&\text{ if }i\le k,\\ g_i&:=z_i&\text{ if }i> k, \end{aligned} \tag{1} \end{equation*} where $$\begin{equation*} s_1:=\sum_1^k z_j^2,\quad s_2:=\sum_{k+1}^n z_j^2. \end{equation*}$$ The Jacobian determinant of the inverse transformation is $$\begin{equation*} J=\det(cM)=c^k\det M, \end{equation*}$$ where $$\begin{equation*} c:=s_2^{1/2}(1-s_1)^{-3/2},\quad M:=(1-s_1)I_k+UU^T, \end{equation*}$$ $$I_k$$ is the $$k\times k$$ identity matrix, and $$U:=[z_1,\dots,z_k]^T$$.

Write $$U=|U|Qe_1$$, where $$|U|=\sqrt{s_1}$$ is the Euclidean norm of $$U$$, $$Q$$ is some orthogonal matrix, and $$e_1:=(1,0,\dots,0)$$. Then it is clear that the matrix $$M$$ is similar to $$N:=(1-s_1)I_k+|U|^2e_1e_1^T=(1-s_1)I_k+s_1e_1e_1^T$$, whence $$\det M=\det N=(1-s_1)^{k-1}$$. So, $$\begin{equation*} J=s_2^{k/2}(1-s_1)^{-k/2-1}. \tag{2} \end{equation*}$$ Also, the joint pdf of $$(G_1,\dots,G_n)$$ is given by $$\begin{equation*} (2\pi)^{-n/2}\exp\Big\{-\frac12\sum_1^n g_j^2\Big\}. \end{equation*}$$ So, in view of (1) and (2), the joint pdf of $$(Z_1,\dots,Z_n)$$ is given by $$\begin{equation*} f_n(z_1,\dots,z_n) =(2\pi)^{-n/2}\exp\Big\{-\frac12\frac{s_2}{1-s_1}\Big\}s_2^{k/2}(1-s_1)^{-k/2-1}. \end{equation*}$$ So, the joint pdf of $$(Z_1,\dots,Z_k)$$ is given by \begin{align*} f_k(z_1,\dots,z_k)&=\int_{\mathbb R^{n-k}}dz_{k+1}\dots dz_n\,f_n(z_1,\dots,z_n) \\ &=(2\pi)^{-n/2}(1-s_1)^{-k/2-1} \\ &\times \int_{\mathbb R^{n-k}}dz_{k+1}\dots dz_n\,s_2^{k/2}\,\exp\Big\{-\frac12\frac{s_2}{1-s_1}\Big\} \\ &\propto(1-s_1)^{(n-k)/2-1}, \end{align*} because $$s_2=\sum_{k+1}^n z_j^2$$. So, we have the desired result.

More specifically, \begin{align*} f_k(z_1,\dots,z_k) =\frac{\Gamma(n/2)}{\pi^{k/2}\Gamma((n-k)/2)}(1-s_1)^{(n-k)/2-1} \end{align*} for $$s_1=\sum_1^k z_j^2\in(0,1)$$.

• The density seems not correct. According to Exercise 1.32 in Aspects of Multivariate Statistical Theory by Robb J. Muirhead, the normalization factor should be $\frac{\Gamma(n/2)}{\color{red}{\pi^{k/2}}\Gamma((n - k)/2)}$. Commented Feb 1, 2023 at 6:09
• @Zhanxiong : Thank you for your comment. The constant factor is now fixed. Commented Feb 1, 2023 at 16:57
• If I followed your calculation, the density I can get is $\frac{\Gamma(n/2)}{(n - k)\Gamma((n - k)/2)\pi^{k/2}}(1 - s_1)^{(n - k)/2 - 1}$, still cannot match the published density. Can you elaborate your calculation? Commented Feb 1, 2023 at 20:29
• By the way, can you check where the argument below is problematic: If let $\sigma(S_{n - 1}(a))$ denote surface area of the $n$-dim sphere, then the joint density of $(X_1, X_2, \ldots, X_n)$ is $1/\sigma(S_{n - 1}(1))$. Partition $(X_1, \ldots, X_n)$ into $(T_1, T_2)$, where $T_1 = (X_1, \ldots, X_k)$, then the marginal density of $T_1$ is $f_{T_1}(t_1) = \int_{t_2: t_1't_1 + t_2't_2 = 1} 1/\sigma(S_{n - 1}(1))dt_2 = \frac{\Gamma(n/2)}{2\pi^{n/2}} \sigma(S_{n - k - 1}(\sqrt{1 - t_1't_1}))$. Commented Feb 1, 2023 at 20:47
• (continued): $=\frac{\Gamma(n/2)}{2\pi^{n/2}}\frac{2\pi^{(n - k)/2}}{\Gamma((n - k)/2)}(1 - t_1't_1)^{(n - k - 1)/2} = \frac{\Gamma(n/2)}{\pi^{k/2}\Gamma((n - k)/2)}(1 - t_1't_1)^{(n - k - 1)/2}$. This is $(1 - t_1't_1)^{-1/2}$ short from the correct answer. But I can't see where the logic went wrong. Commented Feb 1, 2023 at 20:48

$$\newcommand{\R}{\mathbb{R}} \newcommand{\x}{\mathbf{x}} \newcommand{\X}{\mathbf{X}}$$ Here is yet another solution, which is partly informal but I think not hard to completely formalize. Its advantage is a strong and hopefully convincing appeal to geometric intuition.

Again, we are going to assume that $$(X_1,\dots,X_n)$$ is uniformly distributed on the unit sphere $$\mathbb S^{n-1}$$, rather than on $$\sqrt n\,\mathbb S^{n-1}$$. Let $$\X_j:=(X_1,\dots,X_j)$$ and $$\x_j:=(x_1,\dots,x_j)$$ for $$\x=(x_1,\dots,x_n)\in\R^n$$ and $$j=1,\dots,n-1$$. Let $$|\cdot|$$ denotes the Euclidean norm.

The main point is that the probability density $$p_{n-1}(\x_{n-1})$$ of $$\X_{n-1}$$ at a point $$\x_{n-1}\in\R^{n-1}$$ with $$|\x_{n-1}|<1$$ is proportional to the ratio $$r_{n-1}(\x_{n-1}):=vol_{n-1}(dS)/vol_{n-1}(dA)$$, where $$vol_{n-1}$$ is of course the $$(n-1)$$-volume, $$dA$$ is an infinitesimal neighborhood of the point $$\x_{n-1}$$ in $$\R^{n-1}$$ and $$dS$$ is the preimage of $$dA$$ under the projection of the upper hemisphere $$\mathbb S^{n-1}_+:=\{\x\in\mathbb S^{n-1}\colon\x\cdot e_n\ge0\}$$ onto the closed unit ball in $$\R^{n-1}$$; this projection is given by $$\mathbb S^{n-1}_+\ni(\x_{n-1},u)\mapsto\x_{n-1}$$; here $$e_n:=(0,\dots,0,1)$$ and $$\cdot$$ denotes the dot product. But $$$$r_{n-1}(\x_{n-1})=\frac{vol_{n-1}(dS)}{vol_{n-1}(dA)}=\frac1{\cos\phi},$$$$ where $$\phi$$ is the angle between the hyperplane $$\R^{n-1}\times\{0\}$$ of $$\R^n$$ and the tangent hyperplane to $$\mathbb S^{n-1}$$ at the point $$(\x_{n-1},\sqrt{1-|\x_{n-1}|^2})\in\mathbb S^{n-1}_+$$; that is, $$\phi$$ is the angle between the corresponding normal vectors $$e_n$$ and $$(\x_{n-1},\sqrt{1-|\x_{n-1}|^2})$$ of these two hyperplanes.

Thus, $$$$p_{n-1}(\x_{n-1})\propto r_{n-1}(\x_{n-1})\propto\frac1{\cos\phi}=(1-|\x_{n-1}|^2)^{-1/2},$$$$ where $$\propto$$ means an equality up to a constant factor, depending only on $$n,k$$.

Now successively integrating $$p_{n-1}(\x_{n-1})$$ ($$n-1-k$$ times) in $$x_{n-1},\dots,x_{k+1}$$ and each time using the formula $$\begin{equation*} \int_0^{b^{1/2}}(b-u^2)^p du=c_p b^{p+1/2} \end{equation*}$$ for real $$b>0$$ and $$p>-1$$ with $$c_p:=\int_0^1(1-u^2)^pdu\in(0,\infty)$$ (so that $$1/2$$ is added to the exponent $$p$$ after such an integration), we see that the joint pdf of $$\X_k=(X_1,\dots,X_k)$$ is given by $$\begin{equation*} p_k(\x_k)\propto(1-|\x_k|^2)^{(n-k)/2-1}\,1_{|\x_k|<1}, \end{equation*}$$ as desired.

Here is a picture, for $$n=3$$, showing the upper hemisphere $$\mathbb S^{n-1}_+$$; a small neighborhood of a point $$\x_{n-1}$$ in the projection of $$\mathbb S^{n-1}_+$$ to the horizontal plane $$\R^{n-1}\times\{0\}$$ of $$\R^n$$; the preimage of that neighborhood under that projection; and the normal vectors of the horizontal plane and the tangent plane to the sphere -- with $$\phi$$ being the angle between these two normal vectors.

Although great answers have already been provided, the one provided below is perhaps be the shortest possible answer to the question.

Let $$x_1,\ldots,x_k \in \mathbb R$$ such that $$r_k^2:=\sum_{i=1}^k x_i^2 < 1$$, and note that the marginal density of $$(X_1,\ldots,X_k)$$ at $$(x_1,\ldots,x_k)$$ equals $$p_n(r_k^2)$$, where

$$p_n(t) \propto \int_0^\infty\ldots \int_0^\infty dx_{k+1}\ldots dx_n\delta(\sum_{i=1}^n x_i^2-t) \,\forall t \in \mathbb R.$$

The Laplace transform w.r.t $$t$$ is given by $$\begin{split} \hat{p}_n(s) = \int_0^\infty e^{-ts}p_n(t)dt &\propto \int_0^\infty e^{-s\sum_{i=1}^n x_i^2}dx_{k+1}\ldots dx_n \\ &= e^{-s\sum_{i=1}^k x_i^2}\int_0^\infty e^{-s\sum_{i=k+1}^n x_i^2}dx_{k+1}\ldots dx_n\\ &= e^{-sr_k^2}\left(\int_0^\infty e^{-sz^2}dz\right)^{(n-k)} \propto e^{-sr_k^2}s^{-(n-k)/2}. \end{split}$$ Evaluating (e.g via mathematica, etc.) the inverse Laplace transform of the last term at $$t=r_k^2$$, we deduce that $$p_n(r_k^2) \propto (1-t)^{(n-k)/2-1}\delta(t-r_k^2)\bigg\rvert_{t=r_k^2} = (1-r_k^2)^{(n-k)/2-1}$$.

I enjoyed thinking about these answers and this is my attempt to put them into (nonrigorous) geometrical terms. Writing the joint density compositionally as

$$p(\mathbf{x}_k \mid |\mathbf{x}| = 1)p(\mathbf{x}_{n-k} \mid \mathbf{x}_k, |\mathbf{x}| = 1) = p(\mathbf{x} \mid |\mathbf{x}| = 1) \propto 1,$$

we want to solve for the first term on the left. But since our density is proportional to a constant, this is just

$$p(\mathbf{x}_k \mid |\mathbf{x}| = 1) \propto \frac{1}{p(\mathbf{x}_{n-k} \mid \mathbf{x}_k, |\mathbf{x}| = 1)}.$$

Accordingly, instead of performing our calculation by integrating out $$X_{k+1} \dots X_{n}$$, we can think about it in terms of the conditional density for sampling $$\mathbf{X}_{n-k}$$, given $$\mathbf{X}_k$$ and the norm constraint $$|\mathbf{X}| = 1$$, denoted in the denominator above.

I propose a two step procedure. First, draw a point uniformly from within the $$n-k-1$$ dimensional ball with radius $$r_{n-k} = \sqrt{1 - |\mathbf{x}_k|^2}.$$ Each such vector has density proportional to $$\left(1 - |\mathbf{x}_k|^2\right)^{-(n-k-1)/2}.$$

This corresponds to the second and higher integrations in the previous answers, whereas here we directly use the formula for the volume of a ball.

Next, the $$(n-k)$$th coordinate must satisfy $$|\mathbf{x}| = 1$$, which is achieved by any point on the circle with radius $$r_{n-k}$$, a set with measure proportional to $$r_{n-k}$$. Proving this is the first integration in the previous answers.

Putting these two steps together and taking the reciprocal gives

$$p(\mathbf{x}_k \mid |\mathbf{x}| = 1) \propto \left(1 - |\mathbf{x}_k|^2\right)^{(n-k)/2 - 1}.$$