Suppose that $X$ is distributed uniformly in the scaled $n$sphere $\sqrt{n} \mathbf{S}^{n1} \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)^{(nk)/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?
6 Answers
$\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^{n1}$, rather than on $\sqrt n\,\mathbb S^{n1}$.
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_{1t<\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_{n1})
\end{equation*}
for some nonnegative continuous function $g\colon\R^{n1}\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,n1$.
Then
\begin{align*}
\int f\,d\mu_t&=\int_{\R^{n1}}d\x_{n1}\,g(\x_{n1})\int_\R du\, 1_{1t<\x_{n1}^2+u^2\le1} \\
&=\int_{\R^{n1}}d\x_{n1}\,g(\x_{n1})
(1+o(1))t\,(1\x_{n1}^2)^{1/2}\,1_{\x_{n1}<1}.
\end{align*}
Also,
\begin{equation*}
\int d\mu_t=\int_{\R^n}d\x\, 1_{1t<\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^{n1}}d\x_{n1}\,g(\x_{n1})
(1\x_{n1}^2)^{1/2}\,1_{\x_{n1}<1}.
\end{equation*}
Thus, the joint pdf of $\X_{n1}=(X_1,\dots,X_{n1})$ is given by
\begin{equation*}
p_{n1}(\x_{n1})\propto(1\x_{n1}^2)^{1/2}\,1_{\x_{n1}<1}.
\end{equation*}
Now successively integrating $p_{n1}(\x_{n1})$ ($n1k$ times) in $x_{n1},\dots,x_{k+1}$ and each time using the formula
\begin{equation*}
\int_0^{b^{1/2}}(bu^2)^p du=c_p b^{p+1/2}
\end{equation*}
for real $b>0$ and $p>1$ with $c_p:=\int_0^1(1u^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)^{(nk)/21}\,1_{\x_k<1},
\end{equation*}
as desired.

$\begingroup$ This is a wonderfully simple calculation and was illuminating. Thank you for the follow up questions, your persistence, and both answers. $\endgroup$ Commented May 8, 2020 at 20:35

3$\begingroup$ @DrewBrady : Thank you for your kind words. However, I feel grateful to Carlo for teaching me something of value, with patience. $\endgroup$ Commented May 8, 2020 at 21:03

$\begingroup$ Do you by chance know how the ratio$\frac{ \int^{b}_{a} (nx^{2})^{(n2)/2}dx}{\int^{\sqrt{n}}_{\sqrt{n}} (nx^{2})^{(n2)/2}dx} = \frac{ \int^{b}_{a} (1x^{2}/n)^{(n2)/2}dx}{\int^{\sqrt{n}}_{\sqrt{n}} (1x^{2}/n)^{(n2)/2}dx}$ is related to marginal density of the first coordinate? $\endgroup$ Commented Feb 13, 2022 at 15:55

$\begingroup$ 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/annalsofprobability/volume1/… $\endgroup$ Commented Feb 13, 2022 at 15:57

1$\begingroup$ @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,(n1)/2$. The Gaussian approximation to the distribution of $U_1$ easily follows from the representation of $U$ in terms of $n$ iid standard Gaussians. $\endgroup$ 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<n$ elements. The first integration gives $$P_{n1}(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_{n2}(X_3,\ldots X_n)\propto\theta\left(1\sum_{j=3}^n X_j^2\right),$$ the third integration $$P_{n3}(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_{nk+1},\ldots X_n)\propto\left(1\sum_{j=nk+1}^n X_j^2\right)^{(nk)/21}\theta\left(1\sum_{j=nk+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_{n1}$, $$\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^{n1} dr \int d\Omega\, f(r,\phi_1,\phi_2,\ldots\phi_{n1})\frac{2}{A_n}\delta(1r^2)$$ $$\qquad=\frac{1}{A_n}\int d\Omega\, f(r=1,\phi_1,\phi_2,\ldots\phi_{n1}).$$ In the last step I used that $\int_0^\infty r^{n1}\delta(1r^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(1s_2X_1^2),$$ $$\qquad=\int_0^\infty\delta(1s_2q)\frac{dq}{\sqrt q}=(1s_2)^{1/2}\theta(1s_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^2x^2)^p\,dx=c_p a^{1+2p}.$$ }

$\begingroup$ Hm. What does it mean that this distribution "is a delta function"? $\endgroup$ Commented May 7, 2020 at 16:24

$\begingroup$ @ Carlo Beenakker  are you sure you understand what a joint distribution is? $\endgroup$– R WCommented May 7, 2020 at 18:51

$\begingroup$ 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(1x_1^2x_2^2)$ and integrate to obtain $P(\phi)=\int_0^\infty P(r\cos\phi,r\sin\phi) rdr\propto\int_0^\infty \delta(1r^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. $\endgroup$ Commented May 7, 2020 at 21:02

1$\begingroup$ 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 stepbystep integration. Can you write down in detail the integrals for at least the first two steps and how you take them? $\endgroup$ Commented May 7, 2020 at 23:29

1$\begingroup$ I think I got it  thank you. Conversations with physicists are not always easy, but usually fruitful. :) $\endgroup$ 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^{n1}$. 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{1s_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}(1s_1)^{3/2},\quad M:=(1s_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=UQe_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:=(1s_1)I_k+U^2e_1e_1^T=(1s_1)I_k+s_1e_1e_1^T$, whence $\det M=\det N=(1s_1)^{k1}$. So, \begin{equation*} J=s_2^{k/2}(1s_1)^{k/21}. \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}{1s_1}\Big\}s_2^{k/2}(1s_1)^{k/21}. \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^{nk}}dz_{k+1}\dots dz_n\,f_n(z_1,\dots,z_n) \\ &=(2\pi)^{n/2}(1s_1)^{k/21} \\ &\times \int_{\mathbb R^{nk}}dz_{k+1}\dots dz_n\,s_2^{k/2}\,\exp\Big\{\frac12\frac{s_2}{1s_1}\Big\} \\ &\propto(1s_1)^{(nk)/21}, \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((nk)/2)}(1s_1)^{(nk)/21} \end{align*} for $s_1=\sum_1^k z_j^2\in(0,1)$.

$\begingroup$ 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)}$. $\endgroup$ Commented Feb 1, 2023 at 6:09

$\begingroup$ @Zhanxiong : Thank you for your comment. The constant factor is now fixed. $\endgroup$ Commented Feb 1, 2023 at 16:57

$\begingroup$ 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? $\endgroup$ Commented Feb 1, 2023 at 20:29

$\begingroup$ 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}))$. $\endgroup$ Commented Feb 1, 2023 at 20:47

$\begingroup$ (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. $\endgroup$ 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^{n1}$, rather than on $\sqrt n\,\mathbb S^{n1}$. 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,n1$. Let $\cdot$ denotes the Euclidean norm.
The main point is that the probability density $p_{n1}(\x_{n1})$ of $\X_{n1}$ at a point $\x_{n1}\in\R^{n1}$ with $\x_{n1}<1$ is proportional to the ratio $r_{n1}(\x_{n1}):=vol_{n1}(dS)/vol_{n1}(dA)$, where $vol_{n1}$ is of course the $(n1)$volume, $dA$ is an infinitesimal neighborhood of the point $\x_{n1}$ in $\R^{n1}$ and $dS$ is the preimage of $dA$ under the projection of the upper hemisphere $\mathbb S^{n1}_+:=\{\x\in\mathbb S^{n1}\colon\x\cdot e_n\ge0\}$ onto the closed unit ball in $\R^{n1}$; this projection is given by $\mathbb S^{n1}_+\ni(\x_{n1},u)\mapsto\x_{n1}$; here $e_n:=(0,\dots,0,1)$ and $\cdot$ denotes the dot product. But \begin{equation} r_{n1}(\x_{n1})=\frac{vol_{n1}(dS)}{vol_{n1}(dA)}=\frac1{\cos\phi}, \end{equation} where $\phi$ is the angle between the hyperplane $\R^{n1}\times\{0\}$ of $\R^n$ and the tangent hyperplane to $\mathbb S^{n1}$ at the point $(\x_{n1},\sqrt{1\x_{n1}^2})\in\mathbb S^{n1}_+$; that is, $\phi$ is the angle between the corresponding normal vectors $e_n$ and $(\x_{n1},\sqrt{1\x_{n1}^2})$ of these two hyperplanes.
Thus, \begin{equation} p_{n1}(\x_{n1})\propto r_{n1}(\x_{n1})\propto\frac1{\cos\phi}=(1\x_{n1}^2)^{1/2}, \end{equation} where $\propto$ means an equality up to a constant factor, depending only on $n,k$.
Now successively integrating $p_{n1}(\x_{n1})$ ($n1k$ times) in $x_{n1},\dots,x_{k+1}$ and each time using the formula \begin{equation*} \int_0^{b^{1/2}}(bu^2)^p du=c_p b^{p+1/2} \end{equation*} for real $b>0$ and $p>1$ with $c_p:=\int_0^1(1u^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)^{(nk)/21}\,1_{\x_k<1}, \end{equation*} as desired.
Here is a picture, for $n=3$, showing the upper hemisphere $\mathbb S^{n1}_+$; a small neighborhood of a point $\x_{n1}$ in the projection of $\mathbb S^{n1}_+$ to the horizontal plane $\R^{n1}\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^2t) \,\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)^{(nk)} \propto e^{sr_k^2}s^{(nk)/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 (1t)^{(nk)/21}\delta(tr_k^2)\bigg\rvert_{t=r_k^2} = (1r_k^2)^{(nk)/21}$.
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}_{nk} \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}_{nk} \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}_{nk}$, 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 $nk1$ dimensional ball with radius $$r_{nk} = \sqrt{1  \mathbf{x}_k^2}.$$ Each such vector has density proportional to $$\left(1  \mathbf{x}_k^2\right)^{(nk1)/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 $(nk)$th coordinate must satisfy $\mathbf{x} = 1$, which is achieved by any point on the circle with radius $r_{nk}$, a set with measure proportional to $r_{nk}$. 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)^{(nk)/2  1}.$$