Integral of product of Hermite polynomials w.r.t marginal distribution of first two-coordinate of random vector on unit-sphere This question is related to: https://math.stackexchange.com/q/4270522/168758

Let $H_n(x) \in \mathbb R[x]$ be the probabilist's $n$th Hermite polynomial. This an $n$th degree polynomial given by the following equivalent formulae (which ever helps)
$$
\begin{split}
H_n(x)  &= n!\sum_{k=0}^{\lfloor n/2\rfloor}\frac{(-1)^k}{2^kk!(n-2k)!}x^{n-2k}\\
H_n(x) &= \frac{1}{\sqrt{2\pi}}\int_{-\infty}^\infty (x+iy)^n e^{-y^2/2}dy\\
H_n(x) &= e^{-D^2/2}x^n,
\end{split}
$$
where $D^2$ is the second-derivative-w.r.t-$x$ differential operator $\dfrac{d^2}{dx^2}$, and $e^{-D^2/2}$ should be seen as a power series in $D^2$.
Let $d$ be a large positive integer, $a$ and $b$ be fixed vectors on the unit $(d-1)$-dimensional sphere $S_{d-1}$, and $X$ be uniformly distributed on $S_{d-1}$. For fixed nonnegative integers $n$ and $m$, define
$$
s_{n,m} = \mathbb E[H_n(X^\top a)H_m(X^\top b)].
$$
Due to rotational-invarfiance of $X$, it is clear that $s_{n,m}$ is a polynomial in $t:=a^\top b$. Let $c_{n,m,k}$ be the coefficient of $t^k$ in $s_{n,m}$.

Question.
For $k \ge 1$, what is a good Big-O upper-bound for $c_{n,m,k}$ in the limit $d \to \infty$.

 A: To find the dependence of $s_{nm}$ on $t=a\cdot b$, we take $a=(t,\sqrt{1-t^2},0,0,\ldots 0)$, $b=(1,0,0,0,\ldots 0)$, so that
$$s_{nm} = \mathbb E[H_n(X^\top a)H_m(X^\top b)]=\mathbb E[H_n(X_1 t+X_2\sqrt{1-t^2})H_m(X_1)].$$
The marginal distribution $P(X_1,X_2)$ of two elements from a vector that is uniformly distributed on the $d$-dimensional unit sphere is given by (see, for example, this calculation)
$$P(X_1,X_2)=\frac{d-2}{2\pi}(1-X_1^2-X_2^2)^{d/2-2},\;\;X_1^2+X_2^2<1,\;\;d\geq 3.$$
Hence we have for $s_{nm}$ the integral expression
$$s_{nm}=\frac{d-2}{2\pi}\int_{0}^{1}rdr\int_0^{2\pi}d\phi\,  (1-r^2)^{d/2-2}H_n\left(rt\cos\phi+r\sqrt{1-t^2}\sin\phi\right)H_m(r\cos\phi).$$
For large $d$ the Hermite polynomials can be expanded around $r=0$, which gives
$$s_{nm}\approx   \frac{\pi}{d}  2^{\frac{1}{2} (m+n-2)} \left(\frac{4 t}{\Gamma \left(-\frac{m}{2}\right) \Gamma \left(-\frac{n}{2}\right)}-\frac{-2 d+m+n}{\Gamma \left(\frac{1}{2}-\frac{m}{2}\right) \Gamma \left(\frac{1}{2}-\frac{n}{2}\right)}\right),\;\;d\gg 1.$$
A: Disclaimer. This post is just to further simplify @Carlo Beenakker's answer and highlight some potential benefits. It would be a very long comment, so I decided to post it here instead.

With an obvious abuse of notation, let us write $H_n:=H_n(0)$, the $n$th Hermite number. For even $n$, one has
$$
\Gamma(1/2-n/2) = \frac{(-4)^{n/2}(n/2)!\sqrt{\pi}}{n!} = \frac{2^{n/2}2^{n/2}(-1)^{n/2}(n/2)!\sqrt{\pi}}{n!} = \frac{2^{n/2}\sqrt{\pi}}{H_n},
$$
and so we deduce that $\dfrac{\pi2^{(m+n)/2}}{\Gamma(1/2-m/2)\Gamma(1/2-n/2)} = H_nH_m$, and
$$
\begin{split}
\dfrac{\pi2^{(m+n)/2}}{\Gamma(-m/2)\Gamma(-n/2)} &= \dfrac{\pi2^{(m+n)/2}}{\Gamma(1/2-(m+1)/2)\Gamma(1/2-(n+1)/2)}=(1/4)H_{m+1}H_{n+1}
\end{split}
$$
Thus, we get the following instructive formula
$$
s_{nm} \approx \begin{cases}H_nH_m+(1/d)\left(H_nH_{m+1} + H_{n+1}H_m\right),&\mbox{ if }n,m\text{ even},\\
(1/d)H_{n+1}H_{m+1}t,&\mbox{ else,}
\end{cases}
\tag{1}
$$
where we have used the fact that $nH_n = -H_{n+1}$ for every integer $n \ge 0$.

Application
To see the importance of rewriting @Carlo's formula in the form (1), consider the following claim (which settles another question here Approximate the singular values of a certain random dot-product kernel matrix (in the sense of El Karoui, Cheng-Singer, etc.))

Claim. If $g$ is twice continuously-differentiable on $(-1,1)$, then
$$
\mathbb E[g'(X^\top a)g'(X^\top b)]=g'(0)^2+\mathcal O(1/d)+\mathcal O(1/d)t.
$$

It should be noted that the above estimate has been obtained here https://mathoverflow.net/a/405773/78539, under the much more restrictive condition that $g$ is $\mathcal C^2$ on $(-1,1)$ and $\mathcal C^6$ at $0$.
Proof of Claim. Under the hypothesis, $g'$ has a pointwise convergence Hermite expansion (thanks to this post https://mathoverflow.net/a/145235/78539)
$$
g'(x) = \sum_{n \ge 0} b_n(g') H_n(x),\,\forall x \in (-1,1).
$$
In particular, $g'(0) = \sum_{n \ge 0\text{ even }}b_n(g') H_n(0)$.
Here, $b_n(g') := \mathbb E_{z \sim N(0,1)}[g(z)H_n(z)]$ is the $n$th Hermite coefficient of $g$. Recall the important formula
$$
b_n(g') = b_{n+1}(g),\,\forall n \ge 0.
\tag{2}
$$
Now, one has
$$
\begin{split}
\mathbb E[g'(X^\top a)g'(X^\top b)] &= \sum_{n \ge 0}\sum_{m \ge 0} b_n(g')b_m(g')\mathbb E[H_n(X^\top a)H_m(X^\top b)]\\
&= \sum_{n}\sum_{m} b_n(g')b_m(g')s_{n,m}\\
&\overset{(1)}{=} \sum_{n,m \ge 0\text{ even }}b_n(g')H_n(0)b_m(g')H_m(0)\\
&\quad+ \mathcal O(1/d)\sum_{n,m}b_n(g')H_nb_m(g')H_{m+1} +b_n(g')H_{n+1}H_m\\
&\quad+ \mathcal O(1/d)t\sum_{n,m \ge 1\text{ odd }}b_n(g')H_{n+1}b_m(g')H_{m+1}\\
&\overset{(2)}{=} \left(\sum_{n \ge 0}b_n(g')H_n\right)^2\\
&\,+ \mathcal O(1/d)\sum_{n,m}b_n(g')H_nb_{m+1}(g)H_{m+1} +b_{n+1}(g)H_{n+1}H_m\\
&\,+ \mathcal O(1/d)t\sum_{n,m \ge 1\text{ odd }}b_{n+1}(g)H_{n+1}b_{m+1}(g)H_{m+1}\\
&= g'(0)^2 + \mathcal O(1/d)g'(0)(g(0)-b_0(g)) + \mathcal O(1/d)t(g(0)-b_0(g))^2\\
&= g'(0)^2 + \mathcal O(1/d) + \mathcal O(1/d) t,
\end{split}
$$
as claimed.
