Solution (source code)

= Solution

Retain the original PDF's full-line integration limits. For $n\geq2$, write $x^ne^{-x^2}\,dx=-x^{n-1}d(e^{-x^2})/2$. <Integration by parts> gives
$$
I_n=\left[-\frac12x^{n-1}e^{-x^2}\right]_{-\infty}^{\infty}
+\frac{n-1}{2}\int_{-\infty}^{\infty}x^{n-2}e^{-x^2}\,dx.
$$
<Gaussian decay dominates every power>, so the endpoint term vanishes, so
$$
\boxed{I_n=\frac{n-1}{2}I_{n-2}}.
$$
This is the <full-line Gaussian moment recurrence>; its initial values differ from the half-line version.