Integral connection with Hermite and Legendre polynomials Show that $$\int\limits_{-\infty}^{+\infty}x^n e^{-x^2} H_n(tx) dx =\sqrt{\pi} n! P_n(t)$$
Case seems rather complex, I'm completely stuck...
 A: We can also use the following approach. By Rodrigues' formula:
$$H_n(x)=(-1)^ne^{x^2}\frac{d^n}{dx^n}e^{-x^2}$$
Using also
$$\frac{d^n}{dx^n}e^{-x^2}=\frac{1}{2\sqrt\pi}\frac{d^n}{dx^n}\int_{-\infty}^\infty e^{-q^2/4+iqx}dq=\frac{1}{2\sqrt\pi}\int_{-\infty}^\infty e^{-q^2/4+iqx}(iq)^ndq$$
we can present the initial integral in the form
$$I=\int\limits_{-\infty}^{+\infty}x^n e^{-x^2} H_n(tx) dx=\frac{(-1)^n}{2\sqrt\pi}\int_{-\infty}^\infty x^ne^{-x^2(1-t^2)}dx\int_{-\infty}^\infty e^{-q^2/4+iqx}(iq)^ndq$$
Making the substitution $x=\frac{u}{\sqrt{1-t^2}}$ and $q=\frac{v}{\sqrt{1-t^2}}$
$$I=\frac{(-1)^n}{2\sqrt\pi}\frac{1}{1-t^2}\int_{-\infty}^\infty\int_{-\infty}^\infty dudv\,\Big(\frac{iuv}{1-t^2}\Big)^ne^{-u^2}e^{-\frac{v^2}{4(1-t^2)}+\frac{iuv\,t}{1-t^2}}$$
$$=\frac{(-1)^n}{2\sqrt\pi}\frac{1}{1-t^2}\frac{\partial^n }{\partial \color{red}{a}^n}\,\bigg|_{a=t}\int_{-\infty}^\infty\int_{-\infty}^\infty dudv\,e^{-u^2}e^{-\frac{v^2}{4(1-t^2)}+\frac{iuv\,\color{red}{a}}{1-t^2}}$$
Making full square and integrating with respect to $u$ and $v$
$$=\frac{(-1)^n}{2\sqrt\pi}\frac{1}{1-t^2}\frac{\partial^n }{\partial a^n}\,\bigg|_{a=t}\int_{-\infty}^\infty dv\,e^{-\frac{v^2}{4(1-t^2)}}\int_{-\infty}^\infty du\,e^{-\big(u-\frac{iva}{4(1-t^2)}\big)^2-\frac{v^2a^2}{4(1-t^2)^2}}$$
$$=\frac{(-1)^n}{2\sqrt\pi}\frac{1}{1-t^2}\frac{\partial^n }{\partial a^n}\,\bigg|_{a=t}\int_{-\infty}^\infty dv\,e^{-\frac{v^2(1-t^2+a^2)}{4(1-t^2)^2}}\sqrt\pi$$
$$=(-1)^n\sqrt\pi\frac{\partial^n }{\partial a^n}\,\bigg|_{a=t}\frac{1}{\sqrt{1-t^2+a^2}}$$
Making the substitution $a=t-x$
$$I=\sqrt\pi\frac{\partial^n }{\partial x^n}\,\bigg|_{x=0}\frac{1}{\sqrt{1-2xt+x^2}}$$
But we got the generating function for Legendre polynomials:
$$\frac{1}{\sqrt{1-2xt+x^2}}=\sum_{n=0}^\infty P_n(t)\,x^n$$
The result follows.
A: Not an answer, but a possible path to a solution.
A very near miss is (source)
$$\int_0^\infty \mathrm e^{-sx}x^{\nu-1}\operatorname{He}_n(x)\mathrm dx=\Gamma(n+\nu)s^{-n-2}~{}_2F_2\left(\begin{matrix}\frac{-n}{2}~,~\frac{-n+1}{2} \\ \frac{-\nu-n}{2},\frac{-\nu-n+1}{2}\end{matrix}~\bigg |~\frac{-s^2}{2}\right)$$
Where $\operatorname{He}_n$ are the statisticians Hermite polynomials and there appears a generalized hypergeometric function. Unfortunately this identity only holds for odd $n$, but, you can perhaps use the recurrence
$$H_{n+1}(x)=2xH_n(x)-2nH_{n-1}(x)$$
To construct a solution for even $n$.
Of course for your case you can take $s=0$. When $\nu\in\Bbb N$ you can probably use a parity argument to extend the integral to the whole real line. Then perhaps there is some way you can change the ${}_2F_2$ to a ${}_2F_1$, then use the identity
$$P_n(x)={}_2F_1\left(\begin{matrix}-n~,~n+1 \\ 1\end{matrix}~\bigg|~\frac{1-x}{2}\right)$$
Converting the 2F2 to a 2F1 will probably be difficult. functions.wolfram.com has a ton of useful identities listed.
