Suppose $X_1,X_2,\ldots,X_n$ are independent $N(0,\sigma^2)$ random variables.
For $i,j\in \{1,2,\ldots,n\}$, consider $$U=\frac{X_iX_j}{\sum_{i=1}^n X_i^2}$$
Provided $n>1$, we know that $U$ has a Beta distribution when $i=j$ :
$$U=\frac{X_i^2/\sigma^2}{\sum_{i=1}^n X_i^2/\sigma^2} \sim \text{Beta}\left(\frac12,\frac{n-1}2\right) \quad,\,i=1,2,\ldots,n$$
What can we say regarding the distribution of $U$ when $i\ne j$? What are the moments of $U$ in this case?
For $n=2$, if we transform $(X_1,X_2)$ to polar coordinates $(R,\Theta)$, then
$$U=\frac{X_1X_2}{X_1^2+X_2^2}=\frac{R^2\cos\Theta\sin\Theta}{R^2}=\frac{\sin(2\Theta)}{2}$$
Since $\Theta$ is uniformly distributed on $(0,2\pi)$, it seems $\sin(2\Theta)$ has an $\text{Arcsine}(-1,1)$ distribution with pdf
$$f(x)=\frac1{\pi \sqrt{1-x^2}}\mathbf 1_{(-1,1)}(x)$$
So that $U$ has pdf
$$f_U(u)=2 f(2u)=\frac2{\pi\sqrt{1-4u^2}}\mathbf1_{\left(-\frac12,\frac12\right)}(u)$$
If $\boldsymbol X=(X_1,X_2,\ldots,X_n)^T$, we can think of $U$ as the product of $i$th and $j$th components of the vector $\frac{\boldsymbol X}{\lVert \boldsymbol X \rVert}$. And we know that $\frac{\boldsymbol X}{\lVert \boldsymbol X \rVert}$ is uniformly distributed on the surface of a unit sphere. I am not sure if this helps in any way.
Best Answer
Note that we may assume $X_i$'s are standard normal. Then
\begin{align*} \mu_k &:= \mathbf{E}\left[ \biggl( \frac{ X_i X_j }{X_1^2 + \cdots + X_n^2} \biggr)^k \right] \\ &= \frac{1}{(2\pi)^{n/2}} \int_{\mathbb{R}^n} \biggl( \frac{ x_i x_j }{x_1^2 + \cdots + x_n^2} \biggr)^k e^{-(x_1^2 + \cdots + x_n^2)/2} \, \mathrm{d}x_1 \cdots \mathrm{d}x_n \\ &= \frac{1}{(2\pi)^{n/2}} \left( \int_{0}^{\infty} r^{n-1}e^{-r^2/2} \, \mathrm{d}r \right)\left( \int_{\mathbb{S}^{n-1}} \omega_i^k \omega_j^k \, \sigma(\mathrm{d}\omega) \right), \end{align*}
where $\mathbb{S}^{n-1} = \{ \mathrm{x} \in \mathbb{R}^n : \|\mathrm{x}\| = 1\}$ is the unit sphere in $\mathbb{R}^n$ and $\sigma$ is the surface measure on $\mathbb{S}^{n-1}$. Then by mimicking the proof of the beta function identity, such as in this article, it is not hard to prove that
$$ \int_{\mathbb{S}^{n-1}} \omega_1^{\alpha_1} \cdots \omega_n^{\alpha_n} \, \sigma(\mathrm{d}\omega) = \begin{cases} 0, & \text{if some $\alpha_i$ is odd,} \\[0.25em] \dfrac{2\Gamma(\beta_1)\cdots\Gamma(\beta_n)}{\Gamma(\beta_1 + \cdots + \beta_n)}, & \text{if all $\alpha_i$ are even and $\beta_i=\frac{1}{2}(\alpha_i + 1)$}. \end{cases} $$
Using this and simplifying the formula, we get
$$ \mu_k = \begin{cases} 0, & \text{if $k$ is odd,} \\[0.25em] \dfrac{(1 \cdot 3 \cdot 5 \cdots (k-1))^2}{\prod_{l=0}^{k-1} (n+2l)}, & \text{if $k$ is even.} \end{cases} $$
For instance, we get
$$ \mu_2 = \frac{1}{n(n+2)} \qquad\text{and}\qquad \mu_4 = \frac{3^2}{n(n+2)(n+4)(n+6)}.$$