Expectation of a function of two entries of an isotropic unit vector $\mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}\![{w_{1}}^{\!4}\,{w_{2}}^{\!4}]$

Question

Problem:
Given a random isotropic unit vector in $\mathbb{R}^p$ for $p\ge2$, we are trying to compute (preferably exactly, otherwise to upper bound): $$\mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}\!\left[{w_1}^{\!4}\,{w_2}^{\!4}\right]\,.$$

Any help would be greatly appreciated!

---

A related expectation we have already derived:
Using an answer to another question, we were able to compute $\mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}\!\left[{w_1}^{\!2}\,{w_2}^{\!2}\right]$ as follows: \begin{align}1 &=\mathbb{E} \left(\sum_{i=1}^p w_i^2\right)^2 \\[10pt] & = p(p-1)\cdot {\mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}\!\left[{w_{1}}^{\!2}\,{w_{2}}^{\!2}\right]}+p\cdot \mathbb{E} w_1^4 \\[10pt] {\mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}\!\left[{w_{1}}^{\!2}\,{w_{2}}^{\!2}\right]} & =\frac{1-p\cdot \mathbb{E} w_1^4}{p(p-1)}\stackrel{(*)}{=} \frac{1}{p\left(p+2\right)}\,,\end{align} where to show $(*)$ we notice that $w_{1}^{2}=\frac{X}{X+Y}$ where $X\sim\chi^{2}\left(1\right)=\Gamma\left(\frac{1}{2},\frac{1}{2}\right)$ and $Y\sim\chi^{2}\left(p-1\right)=\Gamma\left(\frac{p-1}{2},\frac{1}{2}\right)$, and conclude that $w_{1}^{2}=\frac{X}{X+Y}\sim B\left(\frac{1}{2},\frac{p-1}{2}\right)$ and $\mathbb{E}_{w_1}\left[w_1^4\right] =\operatorname{Var}\left[w_1^2 \right] + \left(\mathbb{E}_{w_1} \left[w_1^2 \right] \right)^2= \dots =\frac{3}{p(p+2)}$.

Answer 1

This is easy to compute using techniques from perturbative quantum field theory, i.e., Wick's Theorem (due to Isserlis) for moments of Gaussian measures.

Consider $$ (2\pi)^{-\frac{p}{2}}\int_{\mathbb{R}^p}d^px\ e^{-\frac{|x|^2}{2}} x_{i_1}\cdots x_{i_n}\ . $$ If $n$ is odd, then it is zero. If $n$ is even, then it is equal to $$ \sum_{M\in P_n^2}\prod_{\{\mu,\nu\}\in M} \delta_{i_{\mu},i_{\nu}}\ , $$ where, following the notations of the above Wikipedia page, $P_n^{2}$ is the set of partitions of $\{1,\ldots,n\}$ into pairs $\{\mu,\nu\}$.

Passing to spherical coordinates and computing an Euler Gamma function, we easily deduce $$ \mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}[w_{i_1}\cdots w_{i_n}] =\frac{1}{p(p+2)\cdots(p+\frac{n}{2}-2)}\sum_{M\in P_n^2}\prod_{\{\mu,\nu\}\in M} \delta_{i_{\mu},i_{\nu}}\ . $$

From this, one easily obtains $$ \mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}[w_1^{2m}w_2^{2m}] =\frac{[1\times 3\times\cdots\times(2m-1)]^2}{p(p+2)\cdots(p+2m-2)} $$ as the special case $i_1=\cdots=i_{2m}=1$, $i_{2m+1}=\cdots=i_{4m}=2$. This is because the ones must be paired with ones, and the twos with twos. We twice get a factor counting perfect matchings of $2m$ elements, which is $1\times 3\times\cdots\times(2m-1)$.

Finally, the answer to the original question is the $m=2$ particular case of the above formula: $$ \mathbb{E}_{\mathbf{w}\sim\mathcal{S}^{p-1}}[w_1^{4}w_2^{4}]= \frac{9}{p(p+2)(p+4)(p+6)} $$ which is very close (and asymtptotic in high dimension) to the upper bound in Alf's comment.

Answer 2

As noted in the previous answer, the joint distribution of $w_1^2$ and $w_2^2$ is the Dirichlet distribution with parameters $1/2,1/2,p/2-1$. Therefore, after some simple calculations we get $$Ew_1^4w_2^4=\frac9{p(p+2)(p+4)(p+6)}. $$